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Localization  of  Deformation  in  Heterogeneous  Solids 
Deformed  at  High  Strain  Rates 


(Abstract) 


Plane  strain  transient  finite  thermomechanical  deformations  of  heat-conducting  particu¬ 
late  composites  comprised  of  circular  tungsten  particulates  in  nickel-iron  matrix  are  analyzed 
using  the  finite  element  method  to  delineate  the  initiation  and  propagation  of  brittle/ductile 
failures  by  the  nodal  release  technique.  Each  constituent  and  composites  are  modeled  as 
strain-hardening,  strain-rate-hardening  and  thermally-softening  microporous  materials.  Val¬ 
ues  of  material  parameters  of  composites  are  derived  by  analyzing  deformations  of  a  rep¬ 
resentative  volume  element  whose  minimum  dimensions  are  determined  through  numerical 
experiments.  These  values  are  found  to  be  independent  of  sizes  and  random  distributions  of 
particulates,  and  are  close  to  those  obtained  from  either  the  rule  of  mixtures  or  microme¬ 
chanics  models. 

Brittle  and  ductile  failures  of  composites  are  first  studied  by  homogenizing  their  material 
properties;  subsequently  their  ductile  failure  is  analyzed  by  considering  the  microstructure. 
It  is  found  that  the  continuously  varying  volume  fraction  of  tungsten  particulates  strongly  in¬ 
fluences  when  and  where  adiabatic  shear  bands  (ASB)  initiate  and  their  paths.  Furthermore, 
an  ASB  initiates  sooner  in  the  composite  than  in  either  one  of  its  constituents. 

We  have  studied  the  initiation  and  propagation  of  a  brittle  crack  in  a  precracked  plate 
deformed  in  plane  strain  tension,  and  a  ductile  crack  in  an  infinitely  long  thin  plate  with 
a  rather  strong  defect  at  its  center  and  deformed  in  shear.  The  crack  may  propagate  from 
the  tungsten-rich  region  to  nickel-iron-rich  region  or  vice-a-versa.  It  is  found  that  at  the 
nominal  strain-rate  of  2000/s  the  brittle  crack  speed  approaches  Rayleighs  wave  speed  in  the 
tungsten-plate,  the  nickel-iron-plate  shatters  after  a  small  extension  of  the  crack,  and  the 
composite  plate  does  not  shatter;  the  minimum  nominal  strain-rate  for  the  nickel-iron-plate 
to  shatter  is  1130/s.  The  ductile  crack  speed  from  tungsten-rich  to  tungsten-poor  regions  is 
nearly  one-tenth  of  that  in  the  two  homogeneous  plates.  The  maximum  speed  of  a  ductile 
crack  in  tungsten  and  nickel-iron  is  found  to  be  about  1.5  km/s. 

Meso  and  multiscale  analyses  have  revealed  that  microstructural  details  strongly  influ¬ 
ence  when  and  where  ASBs  initiate  and  their  paths.  ASB  initiation  criteria  for  particulate 
composites  and  their  homogenized  counterparts  are  different. 
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Chapter  1 
Introduction 


Composites  have  become  common  engineering  materials.  The  prospect  of  producing  “de¬ 
signer”  materials  for  specihc  applications  is  of  great  interest  to  many  engineers.  However, 
the  thermal  and  mechanical  behavior  of  these  composites  is  more  complex  than  that  of 
homogeneous  materials;  many  composites  are  anisotropic  and  exhibit  nonlinear  behavior. 
Experimental  testing  to  determine  the  behavior  of  a  composite  is  tedious  due  to  the  in¬ 
creased  number  of  material  constants  and  the  possibility  of  unexpected  nonlinearities  and 
failure  modes.  Numerical  simulation  provides  the  ability  to  test  a  variety  of  composites  un¬ 
der  various  loading  conditions.  In  this  work,  we  seek  to  better  understand  the  behavior  of 
heterogeneous  materials  at  both  the  micro-  and  the  macro-scales.  The  micro-scale  of  many 
composites  is  far  below  reasonable  length  scales  for  macro-scale  simulations;  accurate  deter¬ 
mination  of  bulk  or  “effective”  properties  that  capture  the  overall  behavior  of  the  composite 
is  highly  desirable.  Furthermore,  failure  modes  and  instabilities  are  not  guaranteed  to  behave 
in  the  same  way  as  in  homogeneous  materials-the  heterogeneity  of  the  material  may  cause 
interesting  effects.  Bridging  the  length  scales  between  the  macro-scale  (at  a  length  scale 
comparable  to  the  material’s  application)  and  the  micro-scale  (at  a  length  scale  comparable 
to  the  microstructure  of  the  material)  allows  one  to  investigate  the  effectiveness  of  the  bulk 
properties  to  capture  heterogeneity-driven  phenomena.  A  greater  understanding  of  these 
phenomena  is  needed  for  critical  uses  of  these  materials. 

In  this  work,  we  focus  on  a  class  of  particulate  composites  known  as  tungsten-heavy  alloys 
(WHAs).  WHAs  are  high  density,  high  strength  materials  that  are  of  interest  in  high  strain- 
rate  loading  scenarios,  most  particularly,  armor  and  armor-penetration  applications.  Due  to 
the  brittleness  of  pure  tungsten,  WHA  was  developed  to  prevent  brittle  fracture/shattering 
under  high  stress;  a  typical  WHA  consists  of  tungsten  particulates  suspended  in  a  nickel-iron 
matrix.  The  particulate  is  typically  spherical  or  ellipsoidal  and  is  randomly  distributed;  thus 
the  overall  composite  is  isotropic.  Unlike  many  typical  composites,  both  constituents  are 
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capable  of  thermoelastoviscoplastic  deformat ions-t he  stiffness  of  tungsten  and  nickel-iron 
are  of  the  same  order  of  magnitude,  but  their  overall  response  when  considering  strain-  and 
strain-rate  hardening,  as  well  as  thermal  softening,  are  decidedly  different. 

The  microstructure  of  a  particulate  composite  lends  itself  to  the  development  of  a  mate¬ 
rial  with  a  spatial  gradation  of  material  properties,  something  frequently  referred  to  in  the 
literature  as  a  functionally  graded  material  (FGM).  FGMs  offer  designers  several  interesting 
possibilities;  most  particularly,  the  ability  to  develop  a  material  with  desirable  properties 
on  two  surfaces  without  the  need  for  a  discrete  interface.  Here,  we  investigate  the  response 
of  tungsten/nickel-iron  FGMs  under  dynamic  loading;  both  in  the  development  of  material 
instabilities  and  brittle  and  ductile  fracture. 

The  dissertation  is  organized  as  follows:  Ghapter  1  describes  previous  work  presented 
by  others  in  the  literature,  and  our  contribution,  if  any,  in  that  area,  as  well  as  the  model 
considered  here.  Ghapter  2  focuses  on  the  development  of  equivalent  thermoelastoviscoplastic 
properties  for  WHAs.  Ghapter  3  shows  the  investigation  into  the  development  of  adiabatic 
shear  bands  in  tungsten/nickel-iron  FGMs.  Ghapter  4  examines  the  brittle  and  the  ductile 
failure  modes  of  tungsten/nickel-iron  FGMs.  Ghapter  5  summarizes  work  done  on  multi-scale 
analyses  that  attempt  to  bridge  the  length  scales  of  WHAs. 

1.1  Literature  Review 

1.1.1  Determination  of  Effective  Properties 

Determination  of  effective  properties  of  a  composite  from  the  properties  of  the  constituents, 
frequently  referred  to  as  homogenization,  has  been  a  topic  of  research  for  several  decades. 
Much  of  the  past  work  has  been  devoted  to  the  determination  of  elastic  constants  of  com¬ 
posites.  Hashin  and  Shtrikman  [1]  established  upper  and  lower  bounds  for  the  shear  and  the 
bulk  modulus  of  a  composite  composed  of  multiple  linear  elastic  phases.  Mori  and  Tanaka  [2] 
developed  a  concept  of  average  stress  to  account  for  the  interaction  among  a  hnite  concen¬ 
tration  of  linear  elastic  ellipsoidal  particles  embedded  in  a  linear  elastic  matrix;  Weng  [3] 
implemented  this  technique  for  isotropic  constituents,  giving  closed  form  expressions  for  the 
shear  and  bulk  modulus  of  the  composite.  Ghristensen  and  Lo  [4]  determined  the  effective 
shear  modulus  of  a  three  phase  composite  by  examining  the  shear  deformation  of  concentric 
spheres.  Equating  the  strain  energy  of  the  composite  sphere  to  the  strain  energy  of  a  ho¬ 
mogeneous  sphere  with  unknown  properties  allows  for  analytical  determination  of  the  shear 
modulus;  however,  this  method  does  not  account  for  interaction  between  particulates.  There 
exist  many  methods  for  determination  of  elastic  parameters  of  various  types  of  materials; 
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these  are  summarized  in  Nemat-Nasser  and  Hori  [5]. 

Estimation  of  the  viscoplastic  parameters  of  a  composite  is  much  more  challenging.  Su- 
quet  [6]  determined  the  upper  and  the  lower  bounds  of  the  yield  stress  of  a  particulate 
composite  composed  of  multiple  incompressible  power  law  materials.  For  a  random  disper¬ 
sion  of  isotropic  spherical  particulates  in  an  isotropic  matrix,  Suquet  derived  an  expression 
for  the  effective  yield  stress  of  the  composite  that  relies  on  the  yield  stress  of  the  constituent 
and  the  effective  shear  modulus  of  the  composite.  Willis  [7]  and  Talbot  and  Willis  [8]  es¬ 
tablished  bounds  for  the  stress-strain  curve  of  a  multi-phase  nonlinear  composite  by  using  a 
variational  technique  similar  to  that  of  Hashin  and  Shtrikman  [1] .  Ponte  Castaeda  [9]  devel¬ 
oped  a  second-order  technique  for  homogenization  of  multi-phase  materials  into  a  nonlinear 
homogeneous  material,  given  prior  knowledge  of  the  stress  concentrations  within  the  com¬ 
posite  and  an  admissible  constitutive  relation.  Many  researchers  have  focused  on  composites 
with  a  viscoplastic  matrix  and  an  elastic  particulate/inclusion  (see  Ponte  Castaneda  [10],  Do- 
hgri  and  Friebel  [11],  Dvorak  and  Srinivas  [12]),  as  this  is  more  common  with  epoxy/carbon 
systems  and  metal  matrix  systems  such  as  aluminum/silicon  carbide  (see  Dai  et  al.  [13]); 
however,  these  techniques  do  not  apply  to  a  system  where  both  the  matrix  and  particulate 
may  experience  large  plastic  deformations. 

1.1.2  Adiabatic  Shear  Banding  in  Functionally  Graded  Materials 

Tresca  [14]  observed  hot  lines  now  called  adiabatic  shear  bands  (ASBs)  during  the  hot 
forging  of  a  platinum  bar;  they  were  subsequently  also  observed  by  Massey  [15].  An  ASB  is 
a  narrow  region,  usually  a  few  micrometers  wide,  of  intense  plastic  deformation  that  forms 
after  softening  of  the  material  due  to  the  heat  generated  by  plastic  working  overcomes  the 
strain-  and  strain-rate  hardening  effects.  Even  though  heat  conduction  plays  a  signihcant 
role  in  determining  the  width  of  an  ASB  and  possibly  spacing  between  adjacent  bands,  the 
adjective  adiabatic  is  used  to  signify  that  they  form  rapidly  in  a  body  deformed  at  a  high 
strain-rate  and  there  may  not  be  enough  time  for  the  heat  to  be  conducted  away.  Their 
study  is  important  because  they  precede  ductile  fracture.  The  research  activity  in  ASBs 
increased  subsequent  to  their  observation  by  Zener  and  Hollomon  [16]  during  the  punching 
of  a  hole  in  a  steel  plate.  They  also  postulated  that  a  material  point  becomes  unstable 
when  thermal  softening  equals  the  combined  hardening  due  to  strain-  and  strain-rate  effects. 
Several  earlier  investigations,  e.g.  see  Recht  [17],  have  used  this  criterion  to  determine  the 
average  strain  at  the  time  of  initiation  of  an  ASB  and  rank  materials  according  to  their 
susceptibility  to  shear  banding;  under  identical  loadings  a  material  in  which  an  ASB  forms 
sooner  is  more  susceptible  to  shear  banding.  Bai  [18]  studied  simple  shearing  deformations 
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of  a  homogeneous  thermoviscoplastic  body,  neglected  elastic  deformations,  perturbed  the 
homogeneous  solution  of  the  governing  equations,  and  hypothesized  that  an  ASB  forms  when 
inhnitesimal  perturbations  grow.  Batra  and  Chen  [19]  have  shown  that  this  approach  gives 
the  same  result  as  that  proposed  by  Zener  and  Hollomon  [16].  Experimental  observations 
of  Marchand  and  Duffy  [20]  on  the  torsion  of  thin-walled  tubes  have  revealed  that  an  ASB 
forms  much  later  than  when  the  shear  stress  peaks  and  is  accompanied  by  a  rapid  drop 
in  the  torque  required  to  deform  the  tube.  Numerical  solutions  of  simple  shearing,  plane 
strain  and  axisymmetric  problems  [21]  and  torsional  deformations  of  a  thin-walled  tube  [22] 
have  conhrmed  that  the  load  carrying  capacity  of  a  member  drops  signihcantly  upon  the 
formation  of  an  ASB.  Much  of  the  literature  on  ASBs  can  be  found  in  two  books  [23,24],  a 
review  paper  [25],  and  special  issues  of  journals  [26-30]. 

One  of  the  earlier  papers  on  FGMs  is  of  Hasselman  and  Youngblood  [31]  who  showed  that 
the  thermal  stress  resistance  of  a  structural  ceramic  can  be  enhanced  by  properly  grading 
the  thermal  conductivity.  Most  works  on  FGMs  have  analyzed  thermoelastic  deformations 
of  plate-like  structures,  e.g,  see  [31-42]  and  those  on  ASBs  have  focused  on  analyzing  their 
initiation  and  propagation  in  a  homogeneous  body.  We  note  that  Batra  and  Zhu  [43-45] 
studied  ASBs  in  a  laminated  composite  body.  Here  we  analyze  the  problem  for  a  FGM 
comprised  of  two  constituents  with  the  compositional  prohle  and  hence  material  properties 
varying  continuously  in  a  cross-section.  A  goal  is  to  find  compositions  that  either  enhance 
or  delay  the  initiation  of  an  ASB.  Whereas  in  a  monolithic  body,  deformations  become 
inhomogeneous  due  to  the  interaction  of  waves  reflected  from  boundaries  with  the  incident 
waves,  in  an  FGM  the  additional  source  of  nonuniform  deformations  is  the  spatial  variation 
of  material  properties.  In  an  FGM  comprised  of  particulates  in  a  matrix,  waves  are  also 
reflected  from  numerous  interfaces  between  the  two  phases.  The  interaction  among  waves 
and  the  sudden  change  in  properties  of  two  constituents  give  inhomogeneous  deformations 
of  the  body.  However,  here  first  a  FG  body  with  continuous  variation  of  material  properties 
is  considered.  Thus  the  wave  speed  varies  continuously  from  point  to  point  and  waves  are 
reflected  only  from  the  boundaries.  The  deformations  are  inherently  inhomogeneous.  It  is 
found  that  the  compositional  prohle  noticeably  affects  the  time  of  initiation,  the  direction 
of  propagation,  and  the  material  point  from  where  an  ASB  forms  relative  to  that  in  a 
homogeneous  body  comprised  of  each  constituent.  Furthermore,  an  ASB  forms  sooner  in 
a  FGM  than  in  either  one  of  the  two  constituents.  Thus  materials  of  desired  susceptibility 
to  adiabatic  shear  banding  can  be  designed.  In  Ghapter  5  we  analyze  thermomechanical 
deformations  of  a  heterogeneous  material  and  compare  its  response  to  that  of  an  equivalent 
homogenized  body. 
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1.1.3  Brittle  and  Ductile  Fracture 

Modeling  crack  propagation  during  the  solution  of  a  transient  problem  by  the  hnite  element 
method  (FEM)  is  very  challenging  since  the  crack  initiation  point  and  its  path  are  to  be 
determined  as  a  part  of  the  solution  of  the  problem.  Three  strategies  often  used  to  ana¬ 
lyze  fracture  are:  (i)  introducing  cohesive  elements  along  inter-element  boundaries  that  are 
weak  in  shear  and  tension  but  very  strong  in  compression,  (ii)  representing  a  crack  as  two 
traction-free  surfaces  by  placing  two  coincident  but  unconnected  nodes  at  the  crack  initiation 
point  and  relieving  tractions  on  the  newly  created  crack  surfaces,  and  (hi)  reducing  elastic 
constants  and  stresses  developed  in  the  failed  region  to  zero  and  virtually  eliminating  these 
elements  from  the  analysis  of  the  problem.  Each  of  these  techniques  has  its  advantages  and 
disadvantages.  Technique  (i)  has  often  been  used  to  simulate  brittle  failure  wherein  strains 
induced  are  small  and  very  little  plastic  deformations  occur.  When  used  to  study  ductile 
failure,  it  does  not  allow  for  the  incorporation  of  frictional  forces  between  two  sliding  surfaces 
formed  during  the  failure  process.  Whereas  technique  (iii)  can  accommodate  large  plastic 
deformations,  it  also  does  not  permit  the  consideration  of  sliding  between  adjacent  surfaces 
surrounding  the  failed  region.  Wang  and  Nakamura  [46]  have  discussed  four  techniques  for 
simulating  material  failure;  however  technique  (ii)  is  not  one  of  them.  They  adopted  tech¬ 
nique  (i)  to  analyze  failure  in  TiB/Ti  EG  ceramic/metal  plate  in  which  brittle  failure  ensues 
at  rather  small  strains.  They  modeled  the  plate  material  as  elastic-plastic  and  considered 
inertia  forces.  An  accurate  modeling  of  the  crack  path  and  hence  its  speed  of  propagation 
requires  a  very  hue  EE  mesh. 

Most  studies  (e.g.,  see  [47-54])  on  fracture  of  FGMs  have  considered  their  static  linear 
elastic  deformations.  Following  an  earlier  study  of  a  propagating  crack  by  Atkinson  and  List 
[55],  there  have  been  at  least  three  studies  [56-58]  on  analyzing  propagating  cracks  in  FGMs. 
Recently,  Lee  [57]  developed  displacement  and  stress  helds  ahead  of  a  crack  propagating  at  a 
uniform  velocity  in  a  FGM  with  shear  modulus  either  varying  linearly  ahead  of  the  crack-tip 
and  Poisson’s  ratio  and  mass  density  kept  constant  or  the  shear  modulus  and  mass  density 
varying  exponentially  but  Poisson’s  ratio  kept  constant.  Works  that  have  analyzed  the 
elastic-plastic  fracture  of  FGMs  include  that  of  Jin  and  Noda  [48]  who  showed  that  the  HRR 
singularity  exists  near  a  stationary  crack-tip  in  a  power-law  hardening  material  provided  the 
yield  stress  is  continuous  and  piecewise  continuously  differentiable  function  of  the  effective 
plastic  strain.  Tvergaard  [58]  modeled  the  interface  between  two  power-law  materials  with  a 
graded  elastic-plastic  layer  and  investigated  crack  growth  in  the  FGM  layer  using  a  cohesive 
zone  model  with  constant  peak  cohesive  traction  and  cohesive  energy  density.  Jin  and 
Dodd  [59]  have  analyzed  the  crack  growth  resistance  behavior  of  a  ceramic/metal  FGM  with 
the  background  material  modeled  by  the  J2  flow  theory  and  undergoing  3-dimensional  quasi- 
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static  plastic  deformations.  They  simulated  crack  propagation  with  a  cohesive  zone  model 
having  six  material  dependent  parameters.  Small  plastic  deformations  are  concentrated  near 
the  crack-tip  and  the  crack  growth  resistance,  calculated  by  the  J-integral,  was  found  to 
increase  signihcantly  with  the  crack  extension  due  to  gradation  in  material  properties.  Wang 
and  Nakamura  [46]  used  the  cohesive  zone  model  to  analyze  transient  inhnitesimal  elastic- 
plastic  deformations  of  a  dynamically  loaded  FGM  plate.  They  simulated  brittle  fracture 
and  found  that  the  energy  of  separation  constituted  a  small  part  of  the  total  fracture  energy. 
The  former  was  computed  from  the  cohesive  law  and  the  latter  equaled  the  sum  of  the 
separation  energy  and  the  energy  dissipated  due  to  plastic  deformations.  Thus,  a  large 
plastic  dissipation  translates  to  higher  crack  growth  resistance.  However,  no  thermal  effects, 
such  as  heat  conduction  and  material  softening  due  to  heating,  were  considered. 

1.2  Formulation  of  the  Problem 

The  study  of  materials  subjected  to  large  deformations  at  high  strain-rates  requires  careful 
formulation  of  the  governing  equations.  The  materials  under  consideration,  tungsten  and 
nickel-iron,  are  both  metallic  materials  that  exhibit  thermo-elasto-viscoplastic  behavior.  The 
occurrence  of  ASBs  during  the  analysis  requires  that  one  consider  large  deformations;  strains 
can  exceed  100%  inside  an  ASB.  Furthermore,  due  to  the  high  strain-rates,  the  constitutive 
theory  employed  must  include  strain-  and  strain-rate  hardening,  as  well  as  thermal  softening. 
Large  plastic  deformations  can  also  give  rise  to  porosity  in  the  material;  the  effect  of  porosity 
on  the  material  parameters  and  the  yield  surface  must  be  considered.  In  the  next  few  sections, 
we  will  outline  the  governing  equations  and  the  FE  formulation  used  to  develop  the  computer 
code  that  is  employed  for  the  studies. 

1.2.1  Governing  Equations 

We  use  rectangular  Cartesian  coordinates  and  the  referential  description  of  motion  to  de¬ 
scribe  transient  coupled  thermomechanical  deformations  of  an  thermo-elasto-viscoplastic 
body  deformed  at  a  high  strain-rate  in  plane  strain  tension.  Deformations  of  each  constituent 
and  the  composite  body  are  governed  by  the  following  equations  expressing,  respectively,  the 
balance  of  mass,  the  balance  of  linear  momentum,  the  balance  of  moment  of  momentum. 
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and  the  balance  of  internal  energy. 


p(l  -  DJ  =  Po(l  -  /o), 

(1.1) 

Po{^  -  fo)vi  =  Tic.,a,  i,j  =  1,2,  a  =  1,2, 

(1.2) 

T  F  =  T  F 

(1.3) 

Po(l  fo)^  Qa,a  “t“ 

(1.4) 

Here  p  is  the  present  mass  density,  /  the  porosity  (i.e.,  the  volume  fraction  of  voids),  J  = 
det  F,  Fia  =  Xi^a  =  dxi/dXa  the  deformation  gradient,  x  the  present  position  at  time  t  of 
a  material  particle  located  at  the  place  X  in  the  reference  conhguration,  T  the  hrst  Piola- 
Kirchhoff  stress  tensor,  e  the  specihc  internal  energy,  Q  the  present  heat  flux  measured  per 
unit  reference  area,  v  the  velocity  of  a  material  particle,  a  superimposed  dot  indicates  the 
material  time  derivative,  and  a  repeated  index  implies  summation  over  the  range  of  the 
index.  Greek  indices  refer  to  coordinates  in  the  reference  conhguration,  and  Latin  indices  to 
coordinates  in  the  present  conhguration. 

We  assume  that  the  strain-rate  tensor  D  dehned  by  Dij  =  -|-  Vj^i)/2,  Vij  =  dvi/dxj, 

has  the  additive  decomposition  into  an  elastic  part  D®,  a  plastic  part  and  a  thermal  part 
a91,  viz.. 


D  =  D^  +  DP  +  d01. 

Equations  (1.1)- (1.5)  are  supplemented  with  the  following  constitutive  relations. 


(1.5) 


<Xij  +  CikWkj  —  (TjkWik  —  T 


1  +  z/ 


(l  +  z/)(l-2z/)' 


e  =  ct6  +  c6  + 


Tia  Jc^iji^F  ')aji 

Qi  ^  ^1  2'^^  )«i) 


(t: 


DP  -  A 


P2P 


2a,. 


3^',  /*/9i/92 


=  ^  -  1  +  2/*/3i  cosh  (  ^  )  -  Plirf  =  0,  al  =  ^abah,  i,j  =  1,  2, 3, 


O':. 


O'!, 


sinh  {  ^  ]  Sij 
2a,  ' 


-y  \^^y 

P  ^  -(<^11  +  *^22  +  <^33)/3,  p  ^  pH{-p  -  0), 


t  P ij  Pi.  P  P^ip 


(1.6) 

(1.7) 

(1.8) 

(1.9) 

(1.10) 

(1.11) 

(1.12) 
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(1  -  f)ayeP 


X  = 


O': 


d(j) 


,  if  0  =  0  and  0  >  0, 


3  when  either  0  <  0  or  0  =  0  and  0  <  0, 

/  =  (1  - /)Z3f,  +  H  ;  H{-p-0), 

S2v27r 


r  = 


/,  /  <  fc, 
r  .  fu  fc 

Jc  “r 


ff-fc 


if  -  fc),  f  >  fc, 


ay  =  {A  +  B{ePX)  [l  +  Cln[^ 


1  - 


6  —  9r 
Orn.  9r 


(1.13) 


(1.14) 

(1.15) 

(1.16) 


The  left-hand  side  of  equation  (1.6)  equals  the  Jaumann  derivative  of  the  Cauchy  stress 
tensor  cr,  Wij  =  {vij  —  Vj^i)/2  is  the  spin  tensor,  E  Young’s  modulus,  z/  Poisson’s  ratio,  a  the 
coefficient  of  thermal  expansion,  6  the  temperature  rise,  6ij  the  Kronecker  delta,  1  the  identity 
tensor,  c  the  specihc  heat,  r  the  thermal  relaxation  time,  and  k  the  thermal  conductivity 
of  the  solid  material.  Batra  and  Jaber  [60]  found  that  the  Jaumann  and  the  Green-Naghdi 
stress  rates  in  equation  (1.6)  give  virtually  identical  results  for  the  times  of  initiation  of  an 
ASB  and  of  brittle  failure  at  points  near  the  surface  of  a  notch  tip  in  a  dynamically  loaded 
prenotched  plate.  This  is  because  elastic  strain-rates  appearing  in  equation  (1.6)  are  small 
as  compared  to  the  plastic  strain-rates  at  points  within  an  ASB.  0  =  0  describes  the  yield 
surface  proposed  by  Gurson  [61]  for  a  porous  material,  p  is  the  hydrostatic  pressure,  and 
f*  the  modihed  value  of  the  porosity  given  by  (1.15).  Gurson’s  yield  surface  is  based  on 
quasistatic  analysis  with  the  matrix  material  modeled  as  rigid  perfectly  plastic  and  obeying 
von  Mises  yield  criterion.  One  should  note  that  in  the  absence  of  porosity,  (1.10)  becomes 
the  classical  von  Mises  yield  surface.  Wang  and  Jiang  [62]  have  considered  inertia  effects. 
Wang  [63]  has  analyzed  transient  deformations  of  a  single  spherical  void  in  a  solid  spherical 
shell  made  of  a  power-law  heat-conducting  viscoplastic  material.  For  a  microporous  material, 
Wang  [63]  derived  an  approximate  expression  for  the  macro-stress  potential  through  an  upper 
bound  approach.  Gonstants  0i  and  02,  introduced  by  Tvergaard  [64],  provide  a  better  £t  of 
results  computed  from  a  FE  analysis  of  the  formation  of  ASBs  in  a  plate  having  an  array  of 
large  cylindrical  voids  with  test  observations,  and  A  is  the  factor  of  proportionality  defined  by 
(1.13);  A  >  0  only  when  the  material  point  is  deforming  plastically,  ay  is  the  current  yield 
stress  of  the  material  whose  dependence  upon  the  effective  plastic  strain  the  effective 
plastic  strain-rate  and  the  temperature  9  is  described  by  the  Johnson-Gook  [65]  relation 
(1.16)  in  which  A,  B,  C,  Jg,  and  m  are  material  parameters,  9r  the  room  temperature  and 
9m  the  melting  temperature  of  the  material.  Parameters  B  and  n  characterize  the  strain 
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hardening  of  the  material,  C  and  the  strain-rate  hardening  and  the  last  factor  on  the  right- 
hand  side  of  (1.16)  its  thermal  softening.  Equation  (1.14)  gives  the  evolution  of  porosity; 
the  hrst  term  on  its  right-hand  side  is  derived  by  assuming  that  the  matrix  is  incompressible 
and  the  elastic  dilatation  is  negligible  as  compared  to  the  plastic  dilatation,  and  the  second 
term  is  the  strain  based  nucleation  of  voids  introduced  by  Chu  and  Needleman  [66].  /2,  S2 
and  En  are  material  parameters;  the  rate  of  nucleation  of  voids  is  highest  when  ef  equals  £n 
and  decays  exponentially  with  the  difference  between  and  H  is  the  Heaviside  step 
function.  Thus  the  second  term  contributes  to  the  evolution  of  porosity  at  a  point  only 
when  the  hydrostatic  pressure  there  is  tensile.  To  account  for  the  coalescence  of  neighboring 
voids,  Tvergaard  and  Needleman  [67]  enhanced  the  porosity,  as  given  by  equation  (1.15), 
after  it  reaches  its  critical  value  fc-  In  equation  (1.15),  //  is  the  porosity  at  ductile  fracture, 
and  fu  =  l/A  is  the  porosity  when  the  yield  surface  has  shrunk  to  a  point.  Equations 
(1.10)  and  (1.16)  imply  that  the  radius  of  the  von  Mises  yield  surface  increases  due  to  strain- 
and  strain-rate  hardening  of  the  material  but  decreases  due  to  the  softening  induced  by 
the  temperature  rise  and  the  increase  in  porosity.  The  degradation  of  material  properties 
due  to  the  damage,  taken  here  synonymous  with  the  porosity,  is  indicated  by  equations  (1.6) 
through  (1.10).  The  affine  variation  with  the  porosity  of  Young’s  modulus,  the  bulk  modulus, 
the  stress-temperature  coefficient,  and  the  heat  capacity  implies  that  the  rule  of  mixtures 
(e.g.  see  Section  2.2.2)  has  been  employed  to  hnd  their  effective  values;  the  expression  for 
the  thermal  conductivity  in  equation  (1.9)  is  due  to  Budiansky  [68].  The  interaction,  if  any, 
among  neighboring  voids  has  been  tacitly  ignored.  Jiang  and  Batra  [69],  among  others,  have 
considered  this  interaction.  The  shrinkage  of  the  yield  surface  due  to  an  increase  in  porosity 
described  by  equation  (1.10)  can  be  seen  by  plotting  the  yield  surface  for  two  different  values 
of  /  while  keeping  other  variables  hxed.  Perzyna  [29]  has  given  a  different  equation  for  the 
evolution  of  porosity. 

For  a  functionally  gradient  or  multi-constituent  body,  all  thermophysical  parameters  may 
vary  with  X. 

Substitution  from  equations  (1.5),  (1.7)  and  (1.9)  into  (1.4)  gives  the  following  hyperbolic 
heat  equation: 

Po(l  -  /o)c(r0  +  9)  =  (1-1’^) 

The  term  JaijD^-  equals  the  heating  due  to  plastic  working  per  unit  volume  in  the  reference 
conhguration;  thus  the  Taylor-Quinney  parameter  has  been  taken  as  1.  Except  for  a  delay  in 
the  time  of  initiation  of  an  ASB  other  results  remain  unaffected  by  a  lower  value  of  the  Taylor- 
Quinney  factor.  The  form  (1.17)  of  the  hyperbolic  heat  equation  is  due  to  Cattaneo  [70]  and 
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Vernotte  [71].  The  thermal  relaxation  time  r  in  it  represents  the  time  required  to  establish 
a  steady  state  of  heat  conduction  in  an  element  suddenly  exposed  to  heat  flux.  According  to 
Chester  [72]  r  equals  “inj pcV^  where  Vq  is  the  speed  of  an  elastic  wave.  Thus  for  a  typical 
steel,  r  =  1  X  10“^^s,  and  r  ~  25  x  10“^^s  for  copper.  Batra  and  Lear  [73],  and  Batra  and 
Chen  [74]  found  that  the  hniteness  of  the  thermal  wave  speed  affects  the  time  of  initiation 
of  an  ASB  in  a  typical  steel  and  the  spacing  between  adjacent  shear  bands  only  when 
r  >  10“®s.  Batra  [75]  considered  higher-order  spatial  and  temporal  gradients  of  temperature 
and  derived  a  heat  equation  that  admits  hnite  speeds  of  thermal  waves.  However,  in  such  a 
material  either  a  thermal  wave  propagates  with  a  hnite  speed  or  the  linearized  problem  has 
a  unique  solution.  Ideally,  one  will  like  to  have  both. 

We  note  that  Batra  and  Kim  [76],  Batra  and  Jaber  [60]  and  Batra  and  Chen  [19]  have 
analyzed  different  aspects  of  shear  banding  with  four  different  thermoviscoplastic  relations, 
namely,  the  Johnson-Cook  [65],  the  Litonski-Batra  (e.g.  see  Batra  [75]),  the  Bodner-Partom 
[77]  and  a  power  law.  These  relations  were  calibrated  to  give  nearly  the  same  effective 
stress  vs.  the  effective  strain  curve  during  homogeneous  deformations  of  the  body.  However, 
during  inhomogeneous  deformations,  each  one  of  the  relations  gave  qualitatively  similar  but 
quantitatively  different  results.  The  decision  to  use  the  Johnson-Cook  relation  here  is  based 
on  the  availability  of  values  of  thermomechanical  parameters  for  different  materials. 


1.2.2  Nondimensionalization  of  Variables 

Let  pij,  di?,  H,  (Jo  and  Or  be  the  reference  mass  density,  the  reference  strain-rate,  the  refer¬ 
ence  length,  the  reference  stress  and  the  reference  temperature  used  to  non-dimensionalize 
quantities.  Then  in  terms  of  nondimensional  variables  indicated  by  the  same  symbols  as 
before,  equations  (1.2)  and  (1.17)  become 


a/(l  -  fo)vi  =  *  =  1,  2,  a  =  1,  2,  (1.18) 

Po(l  -  fo){r9  +  9)  =  -at  (^(^1  -  9,^^  =  1,2,3,  (1.19) 


where 


aj  = 


Pr^rH 


K 


at  — 


9r  = 


O-Q 


(1.20) 


do  pj^cH'^eR  PrCH^' 

aj  and  at  are  nondimensional  measures  of  inertia  and  heat  conduction  effects  respectively. 
For  a  given  material,  inertia  effects  are  directly  proportional  to  the  square  of  the  reference 
strain-rate  and  the  square  of  the  reference  length,  and  heat  conduction  effects  are  inversely 
proportional  to  the  reference  strain-rate  and  the  square  of  the  reference  length. 
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1.2.3  General  Initial  and  Boundary  Conditions 

For  most  analyses,  the  body  is  initially  taken  to  be  at  rest,  with  no  initial  stress  and  at 
uniform  temperature.  There  may  be  an  initial  porosity  distribution  to  simulate  an  initial 
defect.  Thus, 

x(X,  0)  =  X,  v(X,  0)  =  0,  0(X,  0)  =  00,  ^(X,  0)  =  0,  p(X,  0)  =  po, 

^(X,  0)  =  0,  £^(X,  0)  =  0,  /(X,  0)  =  /o(X),  X  G  fi.  (1.21) 

Here  hi  is  the  region  occupied  by  the  body  in  the  reference  configuration. 

Let  the  boundary  dfl  be  divided  into  two  parts,  T^  and  T^,.  On  T^  the  positions  of 
material  particles  are  given  as  functions  of  time,  and  the  tractions  on  T*  are  specihed  as 
functions  of  time.  The  entire  boundary  is  taken  to  be  thermally  insulated.  Thus, 

TiaNa  =  ii  on  Tt,  Xi  =  Xi  on  r„,  QaNa  =  0  on  dVt.  (1-22) 

All  of  the  simulations  conducted  in  Chapters  2  through  5  obey  this  set  of  initial  and  boundary 
conditions.  More  specihc  boundary  conditions  for  each  individual  case,  such  as  the  location 
of  T^  and  r„  and  the  velocity  and  traction  specihcations,  are  given  in  the  respective  chapters. 


1.2.4  Semi-discrete  Formulation  of  the  Problem 

Equations  (1.6),  (1.8)  and  (1.3)  imply  that  the  balance  of  moment  of  momentum  (1.3)  is 
identically  satished.  The  present  mass  density  can  be  computed  from  equation  (1.1)  if 
the  deformation  gradient  and  the  current  value  of  the  porosity  are  known.  However,  here 
the  present  mass  density  at  a  node  is  determined  along  with  other  unknowns.  Thus,  the 
dependent  variables  to  be  solved  for  are  x,  p,  /  and  9  and  the  independent  variables  are  X 
and  t.  Equations  (1.18)  and  (1.19)  are  second-order  coupled  non-linear  hyperbolic  partial 
differential  equations  for  x  and  9.  These  can  not  be  written  explicitly  in  terms  of  x  and 
9  since  T  is  given  by  (1.8)  and  cr  by  (1.6)  which  involves  and  9.  In  order  to  solve  the 
problem  numerically  by  the  EE  method,  we  hrst  derive  its  weak  or  variational  formulation. 

Let  w(X)  be  a  smooth  function  that  is  an  analog  of  virtual  velocity.  We  take  the  inner 
product  of  equation  (1.18)  with  w,  integrate  the  resulting  equation  over  the  region  hi  occupied 
by  the  body  in  the  reference  conhguration,  and  use  the  divergence  theorem  to  arrive  at 


aiil  -  fo)viWidVL  =  - 


WitidV. 


'Tt 


(1.23) 
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Let  '^1,  '^2,  •  •  • ,  t/’n  be  FE  basis  functions  defined  on  We  write 

nodes  nodes 

Vi=  ^  'ljjA(X.)vAi{t),  Wi=  ^  'ljjA(X.)cAi,  i  =  1,2.  (1.24) 

A=1  A=1 

Here  v  is  the  vector  of  velocities  of  nodes,  and  c^j’s  are  constants.  Substituting  from  (1.24) 
into  (1.23)  and  exploiting  the  fact  that  the  resulting  equation  must  hold  for  all  choices  of  c’s 
we  get 


M^=  + 

Mab  =  I  a/(l  -  /o)^A^BdO, 
Jn 

F'm  =  [  i’A,aT,^dn,  Flf  =  [  WMT. 
Jn  Jrt 


(1.25) 


In  order  to  derive  a  weak  form  of  equation  (1.19)  we  hrst  introduce  an  auxiliary  variable 


^  =  0, 


(1.26) 


and  follow  the  same  procedure  as  that  used  to  derive  equation  (1.25)  with  the  following 
result. 


0  =  1 

th|  +  =  f"  +  q, 


(1.27) 


where 


Hab  =  /  Po(l  -  fo)'4)A'4’BdVt, 


rpO  _ 

^A  — 


w  (  1  -  -/  )  9^a'lpA,ad^, 


Qa=  /  t/’A'^fr(crD^)(ir2. 


(1.28) 


Note  that  the  natural  boundary  condition  of  zero  heat  flux  has  been  embedded  in  equation 

(1.27). 

We  solve  equation  (1.16)  for  in  terms  of  ay,  and  6  and  derive  its  weak  form  in  the 
same  way  as  before  except  that  the  divergence  theorem  is  not  used.  Recall  that  df  >  0  only 
when  a  material  point  is  deforming  plastically  as  signihed  by  the  satisfaction  of  equation 
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(1.10);  otherwise  =  0.  Weak  forms  of  equation  (1.6),  (1.14)  and 


x  =  v(X,t)  (1-29) 

are  also  derived.  We  thus  get  coupled  nonlinear  ordinary  differential  equations 

d  =  F,  (1.30) 

where  d  is  the  vector  of  unknowns  and  F  is  the  force  vector  that  depends  upon  time  t  and 
d(t).  The  thirteen  unknowns  at  a  node  are  {xi,X2,Vi,V2,<Jii,o'22,o'i2,<^33,  p}- 

1.2.5  Implementation 

A  computer  code  based  on  the  semidiscrete  formulation  was  developed  using  4-node  quadri¬ 
lateral  elements.  Spatial  integration  is  accomplished  over  each  element  using  2x2  Gaussian 
quadrature.  Spatial  variation  in  material  properties  is  embedded  in  the  integration  by  us¬ 
ing  local  values  at  each  quadrature  point.  Natural  boundary  conditions  are  embedded  in 
equation  (1.30);  suitable  modihcations  to  the  equations  are  made  to  enforce  prescribed  ve¬ 
locity  components  on  T^.  Equation  (1.30)  is  then  integrated  using  the  Livermore  Solver 
for  Ordinary  Differential  Equations  (LSODE);  parameter  ME  in  LSODE  is  set  to  10,  which 
implies  an  Adams-Moulton  method.  The  absolute  and  relative  tolerances  in  LSODE  are  set 
to  10“^.  LSODE  adaptively  adjusts  the  timestep  in  order  to  compute  the  solution  to  within 
the  desired  accuracy.  It  should  be  noted  that  LSODE  is  freely  available  on  the  internet  in 
both  Fortran  and  C-|--|-  variants. 

The  computer  code  was  validated  by  comparing  computed  results  for  several  problems 
with  their  published  analytical  and  numerical  solutions  and  was  also  validated  by  the  method 
of  manufactured  solutions.  In  this  method,  body  forces  and  sonrces  of  internal  energy  den¬ 
sity  are  fonnd  for  any  assumed  deformation  and  temperature  helds  so  that  the  governing 
equations  are  satisfied.  These  helds  are  inpnt  into  the  code.  The  computed  solution  should 
match  the  presumed  analytical  solution  to  a  high  degree  of  accuracy;  see  the  material  fol¬ 
lowing  equation  (20)  in  Batra  and  Liang  [78]  for  an  example.  Furthermore,  the  solution  for 
ASB  initiation  and  propagation  in  a  homogeneous  body  coincided  with  Batra  and  Lear’s  [73] 
solution  obtained  with  triangnlar  elements. 

1.2.6  Simulation  of  Crack  Propagation 

In  order  to  simulate  crack  initiation  and  propagation,  we  assume  that  as  soon  as  a  failure 
criterion  is  met  at  a  node,  say  N,  an  additional  node  N*,  coincident  with  N  bnt  not  connected 
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to  it,  is  added  to  the  FE  mesh.  The  node  N*  is  connected  to  the  node  N**  that  has  the  next 
highest  value  of  the  critical  parameter  that  determines  the  crack  initiation.  The  elements 
are  adjusted  such  that  all  nodes  originally  connected  to  N  on  one  side  of  the  newly  formed 
crack  are  connected  to  N*  instead.  Thus  lines  NN*  and  N*N**  overlap  at  the  instant  of 
the  initiation  of  fracture.  Note  that  no  new  element  is  created;  however,  a  node  is  added, 
the  element  connectivity  is  modihed  and  the  number  of  unknowns  is  increased.  Subsequent 
deformations  of  the  body  will  either  move  N  and  N*  apart  thereby  creating  a  gap  between 
them,  or  tend  to  push  them  together  which  may  be  accompanied  by  relative  sliding  between 
them.  In  the  former  case,  surface  tractions  and  the  normal  component  of  the  heat  flux  are 
assumed  to  be  null  on  the  crack  surfaces.  In  the  latter  case,  the  contacting  surfaces  are 
designated  as  master  and  slave,  and  the  minimum  distance  of  a  node  on  the  slave  surface 
from  the  master  surface  is  determined.  If  interpenetration  occurs,  then  a  Lagrange  multiplier- 
based  conjugate  gradient  method  is  used  to  compute  the  minimum  force  to  be  applied  to  the 
slave  node  and  the  master  line  segment  in  order  to  prevent  interpenetration.  One  returns  to 
the  previous  time  step  and  applies  loads  determined  by  the  contact  algorithm;  this  technique 
is  usually  called  predictor-corrector  and  has  been  described  by  Zywicz  and  Puso  [79],  and 
Carpenter  et  ah  [80].  During  contact,  this  algorithm  is  computationally  expensive,  as  it 
requires  each  timestep  to  be  taken  twice.  However,  the  contact  detection  algorithm  is  fast, 
requiring  little  overhead  during  non-contact  timesteps.  Furthermore,  it  is  independent  of 
which  master  segment  a  node  on  the  slave  surface  penetrates,  gives  the  minimum  normal 
tractions  between  contacting  surfaces  (unlike  less-expensive  penalty  methods),  and  permits 
their  relative  sliding.  Batra  and  Lear  [81]  used  a  penalty-type  method  to  compute  the  normal 
force  to  be  applied  to  contacting  nodes  in  order  to  avoid  interpenetration. 
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Chapter  2 


Equivalent  Thermoelastoviscoplastic 
Parameters  for  Tungsten  Heavy  Alloy 
Particulate  Composites 

2.1  Problem  Statement 

In  this  chapter,  we  seek  to  determine  the  effective  thermoelastoviscoplastic  properties  for 
WHA  particnlate  composites.  These  composites  consist  of  nearly  spherical  tungsten  partic¬ 
ulates,  between  30  and  70p.m  in  diameter,  embedded  in  a  nickel-iron  matrix.  We  utilize  a 
representative  volume  element  (RVE)  indicative  of  the  microstructure  of  the  composite;  this 
RVE  is  subjected  to  a  variety  of  thermomechanical  tests  in  order  to  find  stress-strain  curves. 
These  curves  are  used  to  estimate  strain-  and  strain-rate  hardening  and  thermal  softening 
parameters  for  the  composite.  The  effects  of  RVE  size,  particulate  distribution,  particulate 
diameter,  and  volume  fraction  of  the  particulate  are  considered. 

2.2  Numerical  Tests 

2.2.1  Initial  and  Boundary  Conditions 

Utilizing  the  EE  code  described  in  Section  1.2,  a  square  RVE  with  side  length  H  is  subjected 
to  what  will  be  referred  to  as  a  “constant  area  plane  strain  tension  test”.  The  body  is 
taken  to  be  initially  at  rest,  stress-free,  at  a  uniform  temperature,  with  no  initial  porosity 
(as  indicated  in  Section  1.2.3).  Eurthermore,  porosity  evolution  is  assumed  to  be  zero;  thus, 
the  body  is  always  porosity  free,  giving  a  von  Mises  yield  surface.  The  heat  produced  by 
plastic  working  is  also  assumed  to  be  zero.  The  upper  surface  of  the  body  is  given  a  constant 
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Figure  2.1:  Sketch  of  constant-area  plane-strain  deformation  used  for  homogenization 

velocity  vq  in  the  X2  direction  after  a  short  ramp  time  of  Ifis.  The  lower  surface  of  the  body 
is  held  fixed  in  the  X2  direction  and  is  allowed  to  slide  frictionlessly  in  the  Xi  direction. 
The  lateral  surfaces  are  prescribed  a  velocity  such  that  the  area  of  the  rectangle  remains 
a  constant  as  shown  in  Figure  2.1.  This  choice  of  deformation  is  somewhat  arbitrary; 
however,  it  does  guarantee  that  the  shape  of  the  boundary  remains  well-defined,  despite  the 
inhomogeneous  nature  of  the  RVE.  Also,  the  area-preserving  nature  of  the  RVE  is  consistent 
with  the  incompressibility  assumption  made  with  the  von  Mises  yield  surface,  making  this 
deformation  admissible  for  large  plastic  deformations. 

2.2.2  Rule  of  Mixtures 

According  to  this  rule,  the  value  P  of  a  material  parameter  for  a  mixture  comprised  of 
two  constituents  with  volume  fractions  V/  and  Vf  and  values  Pi  and  P2  of  the  material 
parameter  is  given  by 

P  =  V/Pi  +  R/P2  =  (1  -  R/)Pi  +  R/P2.  (2.1) 

It  gives  exact  values  of  the  mass  density  and  the  heat  capacity  and  often  gives  an  upper 
bound  for  values  of  other  material  parameters  for  the  composite. 

2.2.3  Homogenization  Techniques 

Determination  of  elastic  and  viscoplastic  parameters  from  a  series  of  thermomechanical  de¬ 
formations  of  an  RVE  is  a  non-trivial  task.  In  a  macro-scale  analysis,  the  RVE  itself  is  small 
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enough  to  be  considered  a  material  point;  that  is,  one  could  consider  it  having  uniform  mate¬ 
rial  properties  as  well  as  a  single  stress  and  strain  states.  Thus,  an  effective  stress  and  strain 
tensor  for  the  deformation  of  the  RVE  should  be  dehned.  While  analyzing  deformations  of 
a  RVE  it  has  been  assumed  that  the  two  materials  are  perfectly  bonded  together  at  their 
common  interfaces. 

For  the  constant-area  plane-strain  biaxial  tension/compression  deformation  described  in 
Section  2.2.1,  we  dehne  the  Cauchy  stress  tensor  as  follows:  recalling  that  the  length  in  the 
X3  direction  is  taken  to  be  one,  (T22  is  the  normal  load  divided  by  the  present  length  of  the 
top  and  bottom  edges,  an  is  the  normal  load  divided  by  the  present  length  of  the  lateral 
edges,  ai2  is  the  tangential  load  divided  by  the  present  length  for  all  four  edges  (which  is,  not 
surprisingly,  nearly  zero  for  this  deformation),  and  ass  is  fhe  total  load  in  the  direction 
divided  by  the  current  area  perpendicular  to  the  x^  axis. 

Large  deformations  require  the  dehnition  of  a  suitable  strain  tensor.  The  deformation 
described  in  Section  2.2.1  results  in  the  initial  square  being  deformed  into  a  rectangle.  For 
a  homogeneous  material,  this  would  be  a  homogeneous  deformation;  so,  for  homogenization 
purposes,  we  assume  that  Fn  =  H/a,  F22  =  a/H,  and  Fss  =  1,  where  a  is  the  current  X2 
coordinate  of  the  top  surface  (as  shown  in  Figure  2.1),  and  the  strain  tensor  is  then  dehned 
as: 

£ij  =  —  B-j  )/2  (2.2) 

Here  Bij  =  is  the  left  Cauchy-Green  tensor.  Note  that  the  tensor  e  is  equivalent  to 

the  small-strain  tensor  for  inhnitesimal  deformations.  The  von  Mises  effective  stress  can  be 
computed  as  in  equation  (1.10)2,  and  the  effective  strain  is  taken  to  be  Ee  =  \J ^EijEtj.  Initial 
numerical  analyses  are  conducted  such  that  no  heat  is  generated  during  plastic  working;  that 
is  6^  =  0.  Therefore,  the  thermal  and  mechanical  properties  are  divorced  and  can  be  treated 
separately.  Furthermore,  no  porosity  generation  is  allowed,  reducing  Gurson’s  yield  surface 
(equation  (l.lO)i)  to  the  more  familiar  von  Mises  yield  surface. 

Elastic  Parameters 

Determination  of  Young’s  modulus  and  Poisson’s  ratio  is  relatively  simple  using  the  dehni- 
tions  of  the  stress  and  strain  tensors  in  the  previous  section.  For  metals,  deformations  in  the 
elastic  regime  are  relatively  small  (strains  <  1%),  so  one  may  use  classical  linear  elasticity 
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Figure  2.2:  Comparison  of  yield  stress  estimates  from  the  rule  of  mixtures  and  Suquet’s 
method  using  the  Mori-Tanaka  estimate  for  the  shear  modulus 


relations.  For  plane  strain  deformations  of  isotropic  linear  elastic  materials: 


En  =  -  Z/((T22  +  0-33)], 

(2.3) 

^22  =  -j^[(E22  —  ^{(Eu  +  O'Ss)], 

tj 

(2.4) 

=  0  =  -^[0-33  -  z/(aii  -f  0-22)] 

tj 

(2.5) 

With  the  results  from  a  constant-area  plane  strain  tension  test,  equation  (2.5)  can  be  solved 
directly  for  Poisson’s  ratio,  v.  Using  this  result,  one  can  compute  Young’s  modulus  from 
either  (2.3)  or  (2.4).  The  other  offers  an  excellent  check  on  the  isotropy  of  the  particulate 
arrangement. 

Strain-  and  Strain-Rate  Hardening 

Determination  of  the  quasi-static  yield  stress  A  and  the  strain  hardening  parameters  B  and 
n  in  the  Johnson-Cook  relation  is  difficult  from  an  examination  of  an  effective  stress  versus 
effective  strain  curve.  For  a  given  material,  there  are  several  different  sets  of  parameters  that 
produce  very  close  stress-strain  curves.  The  onset  of  yielding  (corresponding  to  an  effective 
stress  of  A  in  a  quasi-static  test)  is  not  readily  apparent  when  examining  the  stress-strain 
curve  or  its  derivative  because  the  strain-rate  hardening  term  in  the  Johnson-Cook  relation 
prevents  a  noticeable  change  in  the  slope  of  the  curve  at  the  yield  strain.  Figure  2.2  compares 
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values  of  A  computed  from  the  rule  of  mixtures  with  that  obtained  from  Suquet’s  [6]  relation, 
with  the  effective  shear  modulus  obtained  from  the  Mori-Tanaka  [2]  relation.  We  have  chosen 
to  estimate  A  of  the  composite  using  the  rule  of  mixtures.  Conducting  a  plane-strain  tension 
test  (as  described  in  Section  2.2)  at  an  effective  strain-rate  close  to  the  reference  strain  rate 
Eq  then  allows  the  determination  of  the  yield  strain  from  the  estimated  value  of  A;  the 
effective  plastic  strain  ef  is  then  dehned  as  the  effective  strain  minus  the  strain  at  yield. 
Plotting  log((Tjy  —  A)  versus  log(e^)  allows  the  determination  of  B  and  n  from  a  simple 
linear  least  squares  £t.  Through  numerical  tests  it  was  shown  that  the  use  of  a  different 
homogenization  technique  for  determining  A,  such  as  the  method  proposed  by  Suquet  [6], 
produces  quantitatively  different  values  for  B  and  n,  but  the  resulting  stress-strain  curve  is 
effectively  the  same. 

The  effective  strain-rate  hardening  coefficient  C  can  be  determined  by  examining  two 
tests  performed  at  different  effective  strain-rates  on  the  same  RVE.  We  have  taken  the 
reference  strain-rate  Eq  to  be  the  same  for  the  two  constituents  and  the  composite.  Given 
the  stress-strain  curve  for  the  two  different  tests,  the  value  of  C  can  be  calculated  using 
equation  (2.6) 

^ 

^  _ 

—  IniA) 

Here  ayi  and  ay2  are  the  stresses  at  an  arbitrary  strain  during  tests  at  strain-rates  of  £i  and 
£2,  respectively. 

Thermal  Softening 

Thermal  softening  in  the  Johnson-Cook  relation  is  governed  by  6m,  the  melting  temperature 
of  the  material,  and  m,  a  thermal  exponent.  Note  that  as  the  temperature  approaches  dm, 
the  stiffness  of  the  material  approaches  zero.  Here  we  assume  that  the  thermal  exponent  m 
equals  1  for  both  constituents  and  the  composite,  as  it  is  nearly  equal  to  1  for  most  materials. 

To  determine  the  effective  melting  temperature  of  the  composite,  we  conduct  plane  strain 
tension  tests  at  different  uniform  temperatures.  As  in  the  other  tests,  plastic  working  does 
not  contribute  to  temperature  rise;  initial  attempts  at  relaxing  this  assumption  resulted  in 
transient,  non-constant  temperature  fields  in  the  RVE  which  were  dependent  on  the  effective 
strain-rate  of  the  test,  making  it  difficult  to  assign  an  effective  temperature  to  the  RVE. 
Given  the  uniform  temperature,  the  ratio  of  stresses  at  a  given  effective  plastic  strain  for 
two  tests  with  different  temperatures  allows  one  to  simply  determine  the  effective  melting 
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temperature  for  the  composite  from  equation  (2.7) 


9 


m 


^02  +  9i 

(Ty2 


o-j,2 


(2,7) 


Here  and  ay2  are  the  stresses  at  an  arbitrary  strain  during  tests  at  temperatures  of  9i 
and  02,  respectively. 

It  is  important  to  note  that  the  “melting  temperature”  of  the  composite  governs  its 
thermal  softening,  but  has  no  physical  signihcance.  The  actual  melting  temperature  of  the 
composite  would  equal  the  lowest  melting  temperature  among  the  constituents. 


Thermal  Conductivity 

The  effective  thermal  conductivity  of  the  composite  is  determined  by  applying  a  constant 
temperature  differential  across  two  opposite  edges  of  the  RVE.  The  average  heat  flux  on  those 
edges  is  divided  by  the  applied  temperature  gradient,  giving  an  effective  thermal  conductivity. 
Ideally,  one  would  have  a  microstructure  that  gives  an  isotropic  thermal  conductivity;  in 
practice,  this  is  difficult.  The  average  of  the  computed  thermal  conductivity  in  the  Xi  and 
X2  directions  is  taken  to  be  the  effective  thermal  conductivity. 


Coefficient  of  Thermal  Expansion 


The  effective  coefficient  of  thermal  expansion  of  the  composite  is  computed  by  constraining 
the  RVE  on  the  edges  in  the  normal  direction  (as  well  as  the  plane  strain  constraint)  and 
applying  a  uniform  temperature  rise  over  the  entire  RVE.  Given  that  the  resulting  stresses 
are  below  the  yield  stress,  the  effective  coefficient  of  thermal  expansion  is  computed  by 
equation  (2.8): 


^kk 

9i^AT* 


(2.8) 


Here,  K  is  the  effective  bulk  modulus  of  the  material  and  AT  is  the  temperature  rise;  K  is 
computed  from  E  and  v.  Note  that  it  is  assumed  that  the  coefficient  of  thermal  expansion 
is  a  scalar,  tacitly  assuming  an  isotropic  material. 


2.3  Computation  and  Discussion  of  Results 

2.3.1  Tungsten  Heavy  Alloys 

The  particulate  composite  that  we  have  considered  is  composed  of  tungsten  particulates 
dispersed  in  a  nickel-iron  matrix.  When  the  volume  fraction  of  tungsten  exceeds  50%,  the 
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composite  is  frequently  referred  to  as  a  tungsten  heavy  alloy  (WHA);  these  materials  are 
being  considered  as  a  replacement  for  depleted  uranium  in  armor  penetration  applications 
due  to  their  high  density  and  strength.  Pure  tungsten  is  undesirable  in  this  application  due 
to  its  brittleness;  the  addition  of  the  nickel-iron  to  the  matrix  provides  some  ductility  while 
allowing  the  composite  to  retain  some  of  the  strength  and  density  of  tungsten.  Listed  in 
Tables  2.1  and  2.2  are  the  properties  of  the  constituents,  as  given  in  Stevens  and  Batra  [82]. 

For  both  materials,  the  thermal  relaxation  time  was  taken  to  be  r  =  10“®  sec.  The 
thermal  softening  exponent  m  was  taken  to  be  1. 

Note  that  nickel-iron  exhibits  signihcantly  higher  strain-  and  strain-rate  hardening  and 
thermal  softening  than  tungsten.  Numerical  tests  have  shown  that  nickel-iron  is  signihcantly 
less  susceptible  to  shear  band  formation  than  tungsten  [83]. 


Table  2.1:  Elastic  and  thermal  properties  for  W  and  NiFe 


Material 

E 

(GPa) 

V 

P 

(kg/m^) 

K 

(W/m  K) 

a 

(xlO-VK) 

c 

(J/kg  K) 

W 

400 

0.29 

19300 

160 

5.3 

138 

NiFe 

255 

0.29 

9200 

100 

15.0 

382 

Table  2.2:  Viscoplastic  properties  for  W  and  NiFe 


Material 

A 

(MPa) 

B 

(MPa) 

n 

C 

^0 

(1/s) 

n 

(K) 

W 

730 

562 

0.0751 

0.0290 

1.4x10-1^ 

1700 

NiFe 

150 

546 

0.2080 

0.0838 

6.7x10-1^ 

1225 

2.3.2  Verification  of  Homogenization  Methods 

In  an  attempt  to  verify  the  plausibility  and  accuracy  of  the  homogenization  methods  de¬ 
scribed  in  Section  2.2.3,  we  examined  deformations  of  a  pure  tungsten  specimen  of  square 
cross  section  of  length  H  =  0.25  mm.  The  reference  strain-rates  listed  in  Table  2.2  are  far  too 
slow  to  allow  numerical  tests  without  signihcant  strain-rate  effects,  so  the  reference  strain- 
rate  needs  to  be  altered.  If  the  reference  strain-rate  above  is  denoted  by  Sq  and  the  new 
reference  strain-rate  is  denoted  by  ii,  then  the  Johnson-Cook  viscoplastic  relation  (equation 
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1.16,  without  thermal  effects)  gives 


a 


y 


{A  + 

{A  +  B{elT) 
{A  +  B{elT) 


rP 


1  +  Cln{^ 
^0 


gP 

l  +  C'/n( 

Ei  £q 


l  +  C/n  (  ^ 


1  +  C  In 


-C^  In 


Assuming  that  C*  <<  1,  neglecting  the  term  containing  allows  A  and  B  to  be  computed 
at  the  new  reference  strain-rate 


A* 


A 


(^1  +  Cln 


B* 


B 


(^  +  Cln 


(2.9) 


This  modihcation  to  the  material  properties  produces  a  slight  error  (due  to  the  omission  of 
the  term);  however,  this  error  is  of  the  order  of  1  to  2  percent,  far  lower  than  the  accuracy 
level  of  the  material  parameters  themselves.  The  modihed  values  of  A  and  B  are  listed  in 
Table  2.3. 


Table  2.3:  Modified  values  of  A  and  B  for  a  reference  strain-rate  of  5000  1/s 


Material 

A 

(MPa) 

B 

(MPa) 

^0 

(1/s) 

W 

1536.9 

1183.2 

5000 

NiFe 

638.4 

2323.6 

5000 

With  the  modihed  material  parameters,  the  pure  tungsten  RVE  was  subjected  to  tests 
with  nominal  axial  strain-rates  of  2500,  5000,  25000,  and  50000  1/s,  and  additional  50000 
1/s  tests  were  conducted  at  temperature  rises  of  100,  200,  500,  and  1000  K.  Additionally, 
uniform  FE  meshes  of  40  by  40,  60  by  60,  and  80  by  80  elements  were  tested  and  gave 
identical  results  for  the  pure  tungsten  sample.  Using  the  data  from  these  tests,  the  effective 
material  parameters  for  tungsten  were  computed  using  the  homogenization  techniques;  the 
results  are  listed  in  Table  2.4,  along  with  the  percent  error  from  the  input  values. 


Table  2.4:  Values  obtained  from  homogenization  procedure  for  pure  tungsten 


E 

(GPa) 

V 

B 

(MPa) 

n 

C 

0m 

(K) 

K 

(W/mK) 

a 

(lO-^/K) 

homog.  values 

405.94 

0.2907 

1189.4 

0.07813 

0.02877 

1708 

160.0 

5.27 

%  error 

1.49% 

0.24% 

0.52% 

4.03% 

1.03% 

0.47% 

~0%o 

0.57% 
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Note  that  the  quasi-static  yield  stress  A  was  assumed  to  be  1536.9  MPa,  which  was  the 
input  value.  The  accuracy  of  the  homogenization  method  is  excellent,  producing  less  than 
4%  error  for  all  of  the  homogenized  values.  Use  of  different  strain-rate  and  temperature  tests 
resulted  in  slight  differences  in  the  computed  values  of  C  and  9m',  it  was  found  that  using  a 
large  ratio  of  strain-rates  and  a  large  temperature  rise  produced  more  accurate  results. 

2.3.3  Representative  Volume  Element  Size 

The  appropriate  size  for  the  RVE  is  dependent  on  the  microstructure  of  the  material  in 
question.  Typical  WHAs  have  tungsten  particulates  with  diameters  of  ~50  /rm  arranged 
randomly,  such  that  the  bulk  material  is  isotropic.  An  RVE  must  capture  a  sufficient  number 
of  particulates  to  give  an  accurate  volume  fraction  of  material  and  give  a  roughly  isotropic 
response. 

In  order  to  investigate  this  effect,  four  RVEs  were  created  with  side  lengths  of  H  =0.15, 
0.25,  0.35,  and  0.50  mm.  The  EE  mesh  density  was  kept  constant,  with  each  element 
measuring  5  /rm  by  5  /rm.  Each  RVE  consists  of  50  /im  tungsten  particulates  distributed 
randomly  to  give  a  volume  fraction  of  tungsten  equal  to  50%.  Each  of  the  RVEs  was  subjected 
to  a  plane-strain  tension  test  and  then  re-tested  after  a  90  degree  rotation,  to  insure  that 
the  distribution  gave  nearly  isotropic  results.  Shown  in  Figure  2.3  are  all  four  RVEs  and  the 
effective  stress  versus  effective  strain  curves  for  each.  All  four  of  the  curves  agree  pointwise  to 
within  1%.  The  0.25  mm  RVE  was  chosen  for  the  remainder  of  the  analyses,  and  the  effective 
stress  versus  the  effective  strain  curve  obtained  from  it  is  referred  to  as  the  standard  stress 
strain  curve  (the  deviations  in  Figure  2.3  are  normalized  with  respect  to  the  0.25  mm  RVE’s 
result). 

2.3.4  Finite  Element  Mesh  Density 

To  test  the  effects  of  the  mesh  density  on  the  homogenization  results,  a  0.25  mm  RVE 
consisting  of  50%  volume  fraction  tungsten  particulates  of  diameter  50  /rm  was  tested  with 
four  different  uniform  meshes:  40  x  40,  50  x  50,  60  x  60  and  80  x  80  elements.  Figure  2.4 
shows  the  RVE  and  the  deviations  of  the  other  three  meshes  from  the  standard  50  x  50  mesh. 
Note  that  the  mesh  density  does  not  play  a  strong  role  in  influencing  the  homogenization 
results.  The  50  x  50  mesh  was  chosen  for  the  remainder  of  the  analyses. 
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Figure  2.3:  50%  volume  fraction  tungsten  RVEs  with  H  =  (a)  0.15  mm,  (b)  0.25  mm, 
(c)  0.35  mm,  and  (d)  0.50  mm,  and  (e)  the  deviation  in  the  effective  stress 
versus  effective  strain  curve  for  H  =  0.15  mm,  0.35  mm,  and  0.50  mm  from 
that  for  H=0.25  mm 
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(b) 

Figure  2.4:  Finite  element  mesh  density  test:  (a)  0.25  mm  square  50%  volume  fraction 
tungsten  RVE,  (b)  the  deviations  from  the  stress-strain  curve  for  the  40  x 
40,  60  X  60,  and  80  x  80  FE  meshes  from  those  for  the  50  x  50  mesh 
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Figure  2.5:  Deviations  in  the  stress- 
strain  curve  from  the  ordered  particulate 
arrangement  for  three  different  random 
distributions  of  50  /xm  particulate 


Figure  2.6:  Deviations  in  the  stress- 
strain  curve  from  the  50  /um  particulate 
distribution  for  the  30  ^m  and  70  ^um 
distributions 


2.3.5  Particulate  Arrangement 

The  arrangement  of  the  particulates  in  the  RVE  could  possibly  have  an  effect  on  the  values 
of  the  homogenized  parameters.  Most  WHAs  exhibit  isotropic  behavior,  so  the  arrangement 
of  the  particulates  must  be  such  that  the  overall  response  is  independent  of  direction.  To 
examine  this  effect,  four  RVEs  with  50%  volume  fraction  of  50  /xm  tungsten  particulates 
were  tested;  the  first  with  the  particulates  centered  on  the  vertices  of  equilateral  triangles 
and  the  remaining  three  with  different  random  dispersions  of  particulates.  The  deviations 
of  stress-strain  curves  for  these  four  cases  are  shown  in  Figure  2.5,  with  the  reference  curve 
taken  to  be  that  for  the  ordered  arrangement.  No  effect  of  particulate  arrangement  is  seen 
in  the  overall  stress-strain  curve. 


2.3.6  Particulate  Size 

Figure  2.6  shows  the  deviation  in  the  stress-strain  curves  from  the  testing  of  three  RVEs 
with  50%  volume  fraction  of  tungsten  particulates,  each  test  with  a  different  diameter  of 
particulate:  30,  50  and  70  /xm.  The  reference  curve  was  taken  to  be  the  50  /xm  particulate 
case.  One  can  see  no  discernible  effect  of  particulate  size  from  this  analysis. 
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2.3.7  Volume  Fraction  of  Constituents 


Of  primary  interest  to  material  designers  and  engineers  is  the  effect  of  the  volume  fraction 
of  the  particulate  on  the  homogenized  material  parameters.  Results  from  Sections  2.3.3  and 
2.3.4  give  an  indication  of  the  appropriate  RVE  size  and  the  FE  mesh  to  study  the  problem. 
Sections  2.3.5  and  2.3.6  demonstrate  that  the  particulate  size  and  distribution  have  little 
effect  on  the  response  of  the  RVE.  With  these  guidelines,  a  0.25  x  0.25  mm  RVE  meshed 
with  2500  uniform  square  elements  is  studied,  using  50  /rm  diameter  particulates  centered  on 
the  vertices  of  equilateral  triangles.  The  distance  between  the  particulates  is  adjusted  to  give 
volume  fractions  between  50  and  80%.  Less  than  50%  volume  fraction  is  of  little  interest  to 
designers;  greater  than  80%  volume  fraction  is  difficult  to  achieve  with  circular  particulates. 
Each  RVE  was  subjected  to  plane-strain  tension  tests  at  nominal  axial  strain-rates  of  5000 
and  50000  1/s.  Additionally,  the  RVE  was  tested  at  50000  1/s  with  temperature  rises  of  200 
and  500K. 

Elastic  Parameters 

Figure  2.7  shows  Young’s  modulus  as  a  function  of  volume  fraction  of  tungsten.  The  rule 
of  mixtures  and  the  Mori-Tanaka  estimates  are  plotted  as  a  reference.  The  Mori-Tanaka 
estimate,  as  shown  by  Weng  [3],  assumes  randomly  distributed  spherical  particulates;  this 
plane-strain  analysis  assumes  particulates  to  be  cylinders.  Despite  the  non-isotropic  arrange¬ 
ment  (if  one  considers  the  direction),  the  homogenization  estimate  is  quite  reasonable. 
The  Poisson’s  ratio  of  the  composite  matched  the  Poisson’s  ratio  of  the  constituents  within 
1%  for  all  volume  fractions. 

Strain  Hardening  Parameters 

The  variation  of  B  and  n  with  volume  fraction  of  tungsten  is  shown  in  Figures  2.8  and  2.9, 
respectively.  The  homogenization  procedure  predicts  slightly  higher  values  for  B  and  n  than 
the  rule  of  mixtures,  but  the  slope  is  nearly  the  same  as  that  of  the  rule  of  mixtures. 

Strain-Rate  Hardening  Parameters 

The  strain-rate  hardening  coefficient  C  was  determined  for  each  of  the  volume  fractions  at 
effective  plastic  strains  of  5%  through  30%.  Figure  2.10  shows  the  value  of  C  at  5%  effective 
plastic  strain  for  each  of  the  volume  fractions.  The  homogenization  results  are  similar  to  the 
rule  of  mixtures.  It  is  important  to  note  that  the  value  of  C  is  different  over  the  range  of 
strains  studied.  Figure  2.11  shows  C*  as  a  function  of  effective  strain  Ee  for  each  of  the  tested 
volume  fractions. 
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B(GPa) 


Figure  2.7:  Effective  Young’s  modulus  as  a  function  of  volume  fraction  of  tungsten  par¬ 
ticulates 


Figure  2.8:  Strain  hardening  coeffi¬ 
cient  B  as  a  function  of  volume  fraction 
of  tungsten  particulates 
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Figure  2.9:  Strain  hardening  exponent 
n  as  a  function  of  volume  fraction  of 
tungsten  particulates 


Figure  2.10:  Strain-rate  hardening  co¬ 
efficient  C  as  a  function  of  volume  frac¬ 
tion  of  tungsten  particulates 


Figure  2.11:  Strain-rate  hardening  ex¬ 
ponent  C  as  a  function  of  effective  plas¬ 
tic  strain  for  different  volume  fractions 
of  tungsten  particulates 


One  possible  explanation  for  the  variation  in  C  is  the  interaction  between  the  particulates 
and  the  matrix.  Nickel-iron  is  much  softer  than  tungsten,  so  it  constitutes  most  of  the 
initial  deformation.  Also,  it  exhibits  more  strain-rate  hardening;  so  initially,  the  strain-rate 
hardening  value  for  the  composite  is  maximum.  As  the  strain  increases,  more  of  the  tungsten 
is  deformed  plastically,  lowering  the  overall  strain-rate  hardening.  This  effect  is  amplihed  by 
our  choice  of  a  constant  value  for  B  and  n  (dictated  by  the  linear  curve  £t);  realistically,  B 
and  n  are  also  functions  of  strain.  This  variation  with  respect  to  strain  is  dictated  by  the 
microstructure  of  the  material  and  is  difficult  to  quantify.  In  order  to  simplify  the  model, 
B,  n,  and  C  are  taken  to  be  independent  of  strain  for  a  given  volume  fraction,  as  given  by 
Figures  2.8,  2.9,  and  2.10. 


Thermal  Softening 

Figure  2.12  shows  the  melting  temperature  as  a  function  of  volume  fraction.  The  values 
indicated  as  /S.9  =  -I-200A'  are  obtained  by  comparing  the  plane  strain  tension  test  at  the 
reference  temperature  to  a  test  with  a  temperature  rise  of  200K.  The  values  indicated  as 
A9  =  -|-500iF  were  obtained  from  a  test  with  a  temperature  rise  of  500K.  Note  that  there 
is  a  slight  difference  between  the  two  lines.  Figure  2.13  shows  the  variation  in  computed 


28 


Figure  2.12:  Melting  temperature  6m 
as  a  function  of  volume  fraction  of  tung¬ 
sten  particulates 


Figure  2.13:  Melting  temperature  9m 
as  a  function  of  effective  plastic  strain 
for  different  volume  fractions  of  tungsten 
particulates 


melting  temperature  versus  strain.  There  is  a  slight  increase  in  melting  temperature  with 
strain,  possibly  for  similar  reasons  as  the  strain-rate  hardening  case  above. 


Thermal  Conductivity 

Figure  2.14  shows  the  thermal  conductivity  as  a  function  of  volume  fraction.  The  top  and 
bottom  surfaces  of  each  RVE  were  held  at  constant  temperature,  with  the  top  surface  held  at 
-I-IOOK  over  the  bottom  surface.  The  average  heat  flux  in  the  X2  direction  was  calculated  at 
multiple  values  of  X2  =  constant  using  the  local  material  parameters;  this  heat  flux  was  nearly 
constant  throughout  the  body  after  reaching  steady-state  conditions.  The  average  heat 
flux  divided  by  the  overall  temperature  gradient  (lOOK/0.25  mm)  gave  the  values  reported 
in  Figure  2.14.  The  homogenized  values  are  compared  with  the  rule  of  mixtures  and  an 
approach  given  by  Hatta  and  Taya  [84]. 

It  should  be  noted  that  each  RVE  was  tested  at  both  normal  orientation  and  at  a  90° 
rotation.  Eor  the  tungsten/nickel-iron  RVEs,  the  rotation  made  little  difference  in  the  results 
because  of  the  small  difference  in  thermal  conductivity  between  the  metals.  Care  must  be 
taken  when  implementing  this  technique  for  systems  with  great  disparity  between  the  thermal 
conductivities  of  the  constituents;  thermally  isotropic  arrangements  of  such  systems  can  be 
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Figure  2.14:  Thermal  conductivity  n 
as  a  function  of  volume  fraction  of  tung¬ 
sten  particulates 

difficult  to  determine. 

Coefficient  of  Thermal  Expansion 

The  coefficient  of  thermal  expansion  as  a  function  of  volume  fraction  of  tungsten  is  given 
in  Figure  2.15.  The  edges  of  each  RVE  were  held  fixed  as  the  body  was  given  a  uniform 
temperature  rise  of  50K.  The  average  values  of  stress  components  were  calculated  as  indicated 
in  Section  2.2.3,  and  the  homogenized  value  of  the  coefficient  of  thermal  expansion  was 
computed  by  equation  (2.8).  It  should  be  noted  that  the  bulk  modulus  used  in  equation 
(2.8)  is  the  homogenized  value  of  bulk  modulus,  obtained  from  the  values  of  Young’s  modulus 
in  Figure  2.7  and  Poisson’s  ratio  of  0.29.  The  average  stress  (J33  is  slightly  different  from 
that  of  the  in-plane  stresses  an  and  (J22  due  to  the  cylindrical  particulates  required  by  the 
plane-strain  condition;  however,  the  impact  on  the  results  is  negligible.  Note  the  excellent 
agreement  between  the  homogenized  results  and  those  of  the  approach  proposed  by  Hashin 
and  Rosen  [85]. 


Figure  2.15:  Coefficient  of  thermal  ex¬ 
pansion  d  as  a  function  of  volume  frac¬ 
tion  of  tungsten  particulates 


Summary  of  Homogenized  Material  Parameters 

Equations  2.10  through  2.18  give  the  homogenized  material  parameters  as  functions  of  the 
volume  fraction  of  tungsten,  Vf.  A  polynomial  least  squares  £t  was  used  with  sufficient  order 
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to  capture  the  trends.  Note  that  these  relations  are  only  valid  in  the  range  50%  to  80% 
volume  fraction  of  tungsten. 


Young's  modulus  (GPa)  :  E  =  241.76  +  155. 03n/  +  4.53nj,  (2.10) 


Poisson's  ratio  :  u  =  0.29,  (2-11) 

Yield  stress  :  A  =  0.6384  +  0.8985n/,  (2.12) 

Strain  hardening  coefficient  (GPa)  :  B  =  2.394  —  1.152^^;,  (2-13) 

Strain  hardening  exponent  :  n  =  0.2126  —  0.1535n/,  (2-14) 

Strain  rate  hardening  coefficient  :  C  =  0.07563  —  0.04794nj,  (2.15) 

Thermal  conductivity  (W/mK)  :  k  =  76.24  +  81.63n/,  (2-16) 

Goefficent  of  thermal  expansion  (xl0“®/K)  :  a  =  14.11  —  9.138n/,  (2-17) 

Melting  temperature  (°G)  :  9m  =  1214.1  +  447. lu/.  (2.18) 


Listed  in  Table  2.5  are  the  maximum  deviations  of  the  homogenized  values  of  each  thermo¬ 
mechanical  parameter  from  the  rule  of  mixtures  and  corresponding  micromechanical  models. 


Table  2.5:  Maximum  deviation  of  the  homogenized  material  parameters  from  the  rule  of 


mixtures  and  micromechanical  moc 

lels 

Quantity 

Max.  %  deviation  from  the 
Rule  of  Mixtures 

Micromec 

Name 

hanical  model 

Max.  %  deviation 

E 

-2.2% 

Mori-Tanaka 

1.2% 

A 

- 

Suquet 

5.9% 

B 

6.2% 

- 

n 

13.7% 

- 

C 

-10.4% 

- 

- 

9m 

-2.7% 

- 

- 

K 

-10.5% 

Hatta-Taya 

-8.1% 

a 

-7.4% 

Hashin-Rosen 

-3.7% 

2.3.8  Accuracy  of  Homogenization 

In  order  to  test  the  accuracy  of  the  homogenization,  a  homogeneous  material  with  material 
parameters  found  from  the  above  homogenization  techniques  was  subjected  to  a  plane-strain 
tension  test  at  50000  1/s  with  a  uniform  temperature  rise  of  200K.  The  stress-strain  curve 
was  compared  with  a  similar  analysis  performed  on  an  RVE.  Shown  in  Figure  2.16  is  the 
result  for  the  80%  volume  fraction  of  tungsten  case.  The  results  of  a  similar  test  on  a  homo- 
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Figure  2.16:  Effective  stress  versus  effective  strain  curve  for  80%  volume  fraction  tungsten 
RVE  and  its  equivalent  homogeneous  material 


geneous  material  with  the  properties  derived  by  the  rule  of  mixtures  is  given  for  comparison. 
Note  that  the  maximum  point-to-point  error  for  the  homogenized  material  over  the  entire 
range  of  strains  is  1.8%.  For  the  7  different  volume  fractions,  the  maximum  point-to-point 
error  was  2.2%,  with  the  average  point  to  point  error  being  less  than  0.5%.  This  level  of 
accuracy  is  encouraging  considering  that  the  RVE  is  not  truly  isotropic  due  to  the  plane 
strain  assumption  and  inertia  effects  were  neglected  in  the  homogenization  procedure. 


Loading  and  Unloading 

During  a  more  complicated  analysis,  material  points  may  be  loaded,  unloaded,  and  possibly 
reloaded.  To  test  the  accuracy  of  the  homogenization  procedure,  a  RVE  containing  60% 
volume  fraction  of  tungsten  and  its  equivalent  homogeneous  sample  were  subjected  to  the 
plane  strain  tension  test  as  described  in  Section  2.2.3.  However,  the  velocity  on  the  upper 
surface  was  given  a  non-constant  prohle,  given  in  equation  (2.10). 


V2=  { 


on  Xo  =  H. 


(2.19) 


Vf)t,  0  <  t  <  ljus, 

Vo,  1  <  t  <  6jus, 

—vot,  6  <  t  <  8jus, 

—Vo,  t  >  8/is 

Here  the  value  of  Vo  =  12.5m/s,  corresponds  to  a  nominal  axial  strain-rate  of  50000  1/s. 
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Figure  2.17:  Effective  stress  ver¬ 
sus  time  for  plane  strain  tension  un¬ 
load/reload  test  for  60%  volume  fraction 
of  tungsten  RVE  and  its  equivalent  ho¬ 
mogeneous  material 


Figure  2.18:  Driving  force  versus  time 
in  simple  shear  test  for  both  60%  volume 
fraction  of  tungsten  RVE  and  its  equiv¬ 
alent  homogeneous  material 


Shown  in  Figure  2.17  is  a  plot  of  effective  stress  versus  time  for  both  the  homogeneous  sample 
and  the  RVE.  Note  the  excellent  agreement,  even  after  the  body  has  been  allowed  to  unload 
and  then  is  reloaded. 


Simple  Shear 

The  homogenized  material  parameters  are  only  accurate  if  the  response  of  the  homogeneous 
material  compares  well  to  that  of  the  RVE  under  different  types  of  loading.  To  test  the 
homogenization  technique,  a  60%  volume  fraction  tungsten  RVE  and  a  homogeneous  sample 
with  the  homogenized  properties  are  subjected  to  a  simple  shear  test.  The  upper  and  lower 
surfaces  are  restrained  in  the  X2  direction  and  the  upper  surface  is  given  velocity  uq  in  the  xi 
direction  while  the  lower  surface  is  given  the  velocity  —vq  in  the  xi  direction.  Here  uq  was 
chosen  such  that  the  effective  strain-rate  is  approximately  25000  1/s.  Figure  2.18  shows  the 
the  time-histories  of  driving  force  (the  force  in  the  xi  direction  on  the  upper  surface)  for  the 
two  cases.  Note  the  excellent  agreement  between  the  two  sets  of  curves. 
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2.4  Conclusions 


We  have  developed  a  technique  to  hud  the  thermoelastoviscoplastic  properties  of  a  metal- 
metal  particulate  composite.  This  technique  was  used  to  hnd  the  thermoviscoplastic  parame¬ 
ters  for  the  Johnson-Cook  relation  for  a  tungsten  heavy  alloy  composed  of  tungsten  particles 
suspended  in  a  nickel-iron  matrix.  It  was  found  that  the  rule  of  mixtures  under-predicts 
the  strain  hardening  parameters  B  and  n,  with  the  under-prediction  of  n  by  as  much  as 
13%.  Furthermore,  the  strain-rate  hardening  coefficient  C  was  shown  to  vary  with  strain, 
possibly  due  to  the  two  constituents  being  strained  at  different  rates.  The  rule  of  mixtures 
over-predicts  C  by  ~10%.  The  thermal  softening  of  the  homogenized  material,  governed  by 
the  pseudo-melting  temperature  9m  was  found  to  be  slightly  less  than  that  predicted  by  the 
rule  of  mixtures.  A  slight  (<  5%)  variation  in  the  melting  temperature  at  different  strains 
and  temperature  ranges  was  observed. 

The  homogenized  parameters  were  verihed  by  subjecting  a  representative  volume  element 
of  the  microstructure  and  the  corresponding  homogeneous  material  to  a  series  of  tests.  A 
tension/compression  test  where  the  specimen  was  unloaded  and  reloaded,  and  a  simple  shear 
test,  showed  that  the  stress  strain  curves  of  the  homogenized  sample  and  the  RVE  were  nearly 
identical. 
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Chapter  3 

Adiabatic  Shear  Bands  in 

•  *1 

Functionally  Graded  Materials 

3.1  Problem  Statement 

The  initiation  and  propagation  of  ASBs  in  FGMs  is  critical  in  understanding  the  failure  of 
this  new  class  of  materials  when  deformed  at  high  strain-rates.  The  FGM  studied  here  is 
a  tungsten/nickel-iron  alloy;  at  the  micro-scale,  this  material  consists  of  tungsten  particles 
suspended  in  a  nickel-iron  matrix.  Here,  we  study  the  material  at  the  macro-scale,  with  the 
thermomechanical  properties  derived  by  the  rule  of  mixtures  from  the  properties  of  the  two 
constituents.  As  seen  in  Ghapter  2,  the  rule  of  mixtures  gives  macroscopic  properties  that 
are  close  to  those  obtained  by  the  RVE  method.  We  seek  to  hnd  the  effects,  if  any,  that  the 
spatial  gradation  of  material  properties  has  on  the  initiation  and  propagation  of  ASBs. 

3.2  Numerical  Tests 

3.2.1  Initial  and  Boundary  Conditions 

Utilizing  the  plane  strain  FE  code  described  in  Section  1.2,  we  examine  the  thermomechanical 
deformations  of  a  FG  prismatic  body  of  square  cross-section  of  side  length  2H.  We  assume 
that  the  compositional  prohle,  and  thus  the  thermomechanical  deformations,  are  symmetric 
about  the  two  centroidal  axes.  Deformations  of  a  quarter  of  the  cross-section,  represented  by 
the  hrst  quadrant  shown  shaded  in  Figure  3.1,  are  analyzed  using  the  boundary  conditions 
given  in  equation  (3.1).  The  left  (segment  AD)  and  bottom  (segment  AB)  surfaces  contain 

^Material  in  this  Chapter  is  from  Ref.  [83] 
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Vo 


t  t  tj  t  le 


Figure  3.1:  Schematic  sketch  of  the  problem  studied 


the  symmetry  conditions  (no  normal  displacements  or  heat  flux).  The  other  lateral  surface 
(segment  BC)  is  taken  to  be  traction-free.  The  top  surface  CD  is  given  velocity  uq  after  a 
short  ramp  time  of  1/is.  As  described  in  Section  1.2.3,  the  body  is  taken  to  be  initially  at 
rest,  with  uniform  temperature,  and  with  zero  initial  porosity.  Thus 


x(X,  0)  =  X,  v(X,  0)  =  0,  9{X,  0)  =  00,  9{X,  0)  =  0,  p(X,  0)  =  po, 
el{X,0)=0,  /(X,0)  =  0,  {X,,X,)e[0,H]x[0,H], 

T21  =  Til  =  0,  Qi  =  0  on  Xi  =  H, 

T21  =  0,  ui  =  0,  Qi  =  0  on  Xi  =  0, 

Ti2  =  0,  V2  =  0,  Q2  =  0  on  X2  =  0, 

Vot,  0  <t  <  1  /is, 

Vo,  t  >  Ips, 


(3.1) 


Ti2  —  0,  Q2  —  0,  V2  — 


on  Xo  =  H. 
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3.2.2  ASB  initiation  criterion 


Batra  and  Rattazzi  [86]  studied  the  initiation  and  propagation  of  an  ASB  in  a  prenotched 
thick-walled  steel  tube  and  found  that  the  choice  of  the  ASB  initiation  criterion  will  affect  the 
predicted  initiation  time.  They  used  four  different  criteria:  (i)  the  effective  plastic  strain  at  a 
point  equals  0.5,  (ii)  the  effective  plastic  strain  at  a  point  equals  1.0,  (hi)  the  effective  stress 
at  a  point  has  dropped  to  90%  of  its  peak  value  at  that  point,  and  (iv)  the  effective  stress  at 
a  point  has  dropped  to  80%  of  its  maximum  value  at  that  point;  in  each  case  the  material 
point  must  be  deforming  plastically.  Criteria  (iii)  and  (iv)  reflect  Marchand  and  Duffy’s  [20] 
observation  that  the  torque  required  to  twist  thin-walled  tubes  drops  precipitously  upon  the 
initiation  of  an  ASB.  Batra  and  Kim  [87]  scrutinized  ASBs  in  twelve  materials  deformed  in 
simple  shear  and  proposed  criterion  (iv)  which  we  use  in  the  present  work. 


3.3  Computation  and  Discussion  of  Results 

3.3.1  Compositional  Profile 

Results  have  been  computed  for  tungsten  (W)  particles  interspread  in  NiFe  matrix  with  the 
volume  fraction,  Vf^NiPei  of  NiFe  given  by  one  of  the  following  two  expressions. 


^f,NiFe 


Type-I 

r  r  <  H 

1  C/,  r  >  H 


^f,NiFe 


Type-II 

C/(l-i).  r<H 
0,  r  >  H. 


(3.2) 


Here  r  is  the  radial  distance  from  the  specimen  centroid,  and  2H  =  10  mm  equals  the  length 
of  a  side  of  the  square  cross-section.  In  both  type-I  and  type-II  FGMs,  material  properties 
in  the  region  r  >  H  are  constants.  In  type-I  FGMs,  the  volume  fraction  of  NiFe  varies  from 
zero  at  the  specimen  centroid  to  cj  at  r  =  77,  and  in  type-II  FGMs,  it  varies  from  cj  at  the 
specimen  centroid  to  zero  at  r  =  77.  The  material  is  homogeneous  in  the  region  r  >  77.  The 
variation  of  the  mass  density  in  these  two  types  of  FGMs  is  depicted  in  Figure  3.2. 

The  primary  reason  for  studying  compositional  prohle  (3.2)  is  that  it  does  not  influence 
the  direction  of  propagation  of  an  ASB  originating  from  the  centroid  of  the  cross-section.  For 
the  same  overall  volume  fraction  of  the  constituents,  the  ASB  initiation  time  may  depend 
upon  the  distribution  of  the  two  constituents.  If  an  ASB  does  not  initiate  from  the  centroid, 
then  the  distribution  (3.2)  of  NiFe  may  affect  the  direction  of  propagation  of  the  ASB. 
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3.3.2  Finite  Element  Mesh 

In  order  to  ascertain  the  effect  of  the  FE  mesh  on  the  ASB  initiation  time,  we  analyzed 
deformations  of  type-I  FGM  with  c/  =  0.3  in  equation  (3.2)  with  three  different  uniform 
FE  meshes.  As  is  clear  from  the  results  summarized  in  Table  3.1,  the  ASB  initiation  time 
decreased  by  2.1%  in  going  from  a  40  x  40  to  a  120  x  120  uniform  FE  mesh,  and  the  CPU 
time  increased  by  a  factor  of  21.  Results  presented  in  subsequent  sections  are  with  a  40  x  40 
uniform  EE  mesh  of  4-node  isoparametric  quadrilateral  elements. 


Table  3.1:  Effect  of  FE  mesh  on  the  ASB  initiation  time  for  a  type-I  EG  material  with 

%  =  0-3- _ 


Uniform 
FE  mesh 

ASB  initiation 
time  (ps) 

CPU  time 
(secs) 

%  change  in 

ASB  initiation  time 

40  X  40 

65.9 

1,133 

— 

80  X  80 

64.8 

6,908 

1.67 

120  X  120 

64.5 

29,242 

2.12 

3.3.3  ASBs  in  Five  Homogeneous  Materials 

For  homogeneous  materials,  Batra  and  Lear  [45]  assumed  the  presence  of  2.5%  or  5%  initial 
porosity  at  the  specimen  centroid  which  decayed  exponentially  to  zero  with  the  distance 
from  the  centroid.  Because  of  the  maximum  initial  porosity,  the  centroid  acts  as  a  nucleation 
site  for  an  ASB.  However,  for  a  composite,  deformations  naturally  become  inhomogeneous 
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because  of  the  spatial  variation  of  material  properties  and  it  determines  where  and  when  an 
ASB  initiates.  For  reference  and  to  delineate  the  effect  of  the  initial  porosity  distribution,  we 
have  listed  in  Table  3.2  values  of  material  parameters,  and  in  Table  3.3  the  nondimensional 
ASB  initiation  times  (or  equivalently  the  average  axial  strain)  for  hve  materials  deformed 
at  an  average  axial  strain-rate  of  5,000/s.  For  H  =  5  mm,  the  axial  speed  Uq  was  increased 
from  zero  to  25  m/s  in  1  /is  and  then  kept  constant.  It  is  clear  that  an  ASB  initiates 
considerably  sooner  in  initially  porous  materials  than  in  those  with  zero  initial  porosity. 
Nondimensional  inertia  factor  aj  is  higher  for  NiFe,  Armco  iron  and  OFHC  copper  than 
that  for  4340  steel  and  tungsten  even  though  the  mass  density  of  tungsten  is  almost  twice 
that  of  steel.  Heat  conduction  effects  are  highest  in  OFHC  copper  than  in  any  of  the  other 
four  materials  studied;  it  is  reflected  in  the  delayed  initiation  of  an  ASB  in  copper.  The 
rather  high  nondimensional  thermal  conductivities  of  NiFe  and  Armco  iron  also  delay  the 
onset  of  ASBs  in  them.  The  nondimensional  thermal  conducitivites  of  tungsten  and  NiFe 
are  nearly  the  same  but  for  tungsten  is  about  l/5th  of  that  for  NiFe.  A  clear  correlation 
among  the  ASB  initiation  time,  aj  and  at  has  not  been  established  yet.  We  note  that  the 
ASB  initiation  time  also  depends  upon  the  nominal  axial  strain-rate. 


Table  3.2:  Values  of  thermomechanical  parameters  for  the  five  materials  studied 


Material 

A 

(MPa) 

B 

(MPa) 

C 

m 

n 

dm 

(K) 

eo 

(1/s) 

4340  steel 

792.19 

509.51 

0.014 

1.03 

0.26 

1793 

10-^ 

Armco  iron 

175.12 

379.99 

0.060 

0.55 

0.32 

1811 

10“^ 

OFHC  Copper 

89.63 

261.64 

0.031 

1.09 

0.31 

1356 

10“^ 

Tungsten 

1536.9 

1183.2 

0.029 

1.00 

0.0751 

1700 

5000 

Nickel-Iron 

638.4 

2328.6 

0.0838 

1.0 

0.208 

1225 

5000 

Material 

a 

(10-Vk) 

K 

(W/mK) 

c 

(J/kgK) 

P 

(kg/m^) 

V 

E 

(GPa) 

* 

aj 

(10-3) 

* 

at 

(10-6) 

4340  steel 

12.3 

38 

477 

7840 

0.27 

210 

6.19 

81.3 

Armco  iron 

12.1 

73 

452 

7890 

0.29 

196 

28.16 

163.8 

OFHC  copper 

17.0 

386 

383 

8960 

0.33 

120 

62.48 

899.9 

Tungsten 

5.3 

160 

138 

19300 

0.29 

400 

7.849 

480.6 

Nickel-Iron 

15.0 

100 

382 

9200 

0.29 

255 

9.007 

227.6 

*Values  for  an  average  axial  strain-rate  of  5,000/s,  H  =  5mm,  and  (Tq  =  A. 


Figure  3.3  exhibits  the  time-history  of  the  applied  axial  load  for  the  hve  materials;  the 
axial  strain  plotted  along  the  abcissa  equals  the  time  multiplied  by  the  nominal  axial  strain- 
rate.  In  each  case  the  initiation  of  an  ASB  is  accompanied  by  a  rapid  drop  in  the  axial 
load;  the  rate  of  drop  of  the  axial  load  is  highest  for  tungsten  than  that  for  any  of  the  other 
four  materials.  In  each  case  the  ASB  initiated  from  the  centroid  of  the  cross-section.  The 
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Table  3.3:  Axial  strain  at  ASB  initiation  in  five  homogeneous  materials  deformed  in  plane 


strain  tension  at  a  nominal  strain-rate  of  5,000/s. 


Material 

Axial  strain  at  ASB  initiation 

Temperature 
rise  (K)  at  ASB 
initiation* 

w/  initial  porosity 

w/o  initial  porosity* 

Tungsten 

0.137 

0.386 

391 

4340  steel 

0.315 

0.615 

315 

Nickel  iron 

- 

0.944 

354 

OFHG  copper 

0.391 

1.045 

157 

Armco  iron 

0.433 

1.095 

253 

*  denotes  present  work 


porosity  at  the  instant  of  the  ASB  initiation  equaled  0.041,  0.087,  0.081,  0.126  and  0.155 
respectively  in  W,  NiFe,  4340  steel,  Fe  and  Cu;  the  temperature  rise  is  listed  in  Table  3.3. 
An  approximate  value  of  the  temperature  rise  is  given  by  ((A  -|-  0.5il(e^)"')e^)/(pc)  since 
the  effect  of  heat  conduction  till  an  ASB  initiates  is  small.  It  also  provides  a  check  on  the 
accuracy  of  computed  results.  For  W,  the  temperature  rise  computed  from  this  relation  at 
the  time  of  initiation  of  an  ASB  equals  368  K. 

3.3.4  ASBs  in  FGMs 

Numerical  experiments  on  type-I  and  type-II  FGMs  yielded  interesting  results  summarized 
in  Table  3.4.  In  type-I  FGMs,  there  is  100%  W  at  r  =  0,  and  the  volume  fraction  of  NiFe 
equals  0  at  r  =  0  and  increases  linearly  with  r  reaching  its  maximum  value  at  r  =  77;  the 
maximum  volume  fraction  of  NiFe  at  r  =  77  varies  from  0.1  to  0.5.  In  each  case,  an  ASB 
initiated  at  the  specimen  centroid  and  propagated  along  a  line  inclined  at  approximately  45° 
to  the  horizontal  axis  in  the  present  conhguration.  The  axial  strain  at  the  instant  of  ASB 
initiation  decreased  when  the  maximum  volume  fraction  of  NiFe  at  r  =  77  was  increased 
from  0.1  to  0.4,  and  then  increased  a  little  when  c/  in  equation  (3.3)  was  increased  from  0.4 
to  0.5.  However,  in  each  case,  the  average  axial  strain  at  the  instant  of  the  ASB  initiation  is 
less  than  that  for  the  pure  W.  It  leads  one  to  conjecture  that  the  addition  to  W  of  a  material 
less  susceptible  to  shear  banding  than  W  increases  tungsten’s  susceptibility  to  shear  banding. 
In  order  to  test  this  conjecture,  numerical  experiments  were  performed  on  type-I  FGMs  with 
NiFe  replaced  by  4340  steel,  Gu  and  Fe.  The  results  of  these  experiments  are  also  included  in 
Table  3.4.  The  spatial  distribution  of  the  volume  fraction  of  W  was  kept  the  same.  Whereas 
in  each  FGM  the  ASB  initiated  sooner  than  that  in  pure  W,  the  decrease  in  the  time  of 
initiation  of  an  ASB  cannot  be  simply  related  to  the  times  of  initiation  of  ASBs  in  the 
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Figure  3.3:  Axial  load  versus  axial  strain  for  five  homogeneous  materials  deformed  in 
plane  strain  tension 


constituents. 

Contours  of  the  effective  plastic  strain  in  type-I  FGMs  are  depicted  in  Figure  3.4.  The 
effective  plastic  strain  within  the  ASB  exceeds  1  with  a  peak  value  of  1.5  at  some  points. 
The  necking  of  the  region  below  the  ASB  is  evident  in  these  plots. 

For  type-II  FGMs  with  the  maximum  volume  fraction  of  NiFe  at  the  centroid,  an  ASB 
initiated  from  a  point  on  the  top  surface  and  propagated  inwards  along  the  two  directions 
of  the  maximum  shear  stress.  These  seem  to  get  reflected  from  the  surfaces  Xi  =  0,H. 
As  shown  in  Figure  3.5,  there  are  several  narrow  regions  of  large  plastic  deformation.  For 
the  hrst  three  FGMs,  the  point  where  an  ASB  initiates  on  the  top  surface  varies  with  the 
compositional  prohle.  For  the  other  two  FGMs,  an  ASB  initiates  from  a  point  within  the 
cross-section.  ASBs  exhibited  in  Figure  3.5  reveal  that  the  eventual  distribution  of  regions  of 
intense  plastic  deformation  is  independent  of  c/.  The  time  of  initiation  of  an  ASB  decreases 
monotonically  with  an  increase  in  the  maximum  volume  fraction,  Cf,  of  NiFe.  However,  the 
rate  of  decrease  of  the  time  of  initiation  decreases  rapidly  as  c/  is  increased  from  0.1  to  0.5. 
In  contrast  to  FGMs  of  type-I,  necking  occurs  in  the  upper  part  rather  than  in  the  lower 
portion  of  the  cross-section. 
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Table  3.4:  Axial  strain  and  point  of  initiation  of  ASBs  in  FGMs 


Specimen  type 
and  constituents 

Max.  uj  of 
NiFe 

Axial  strain  at 
ASB  initiation 

Ref.  Location  of  ASB  initiation 

X 

Y 

W/NiFe  Type-I 

0.1 

0.371 

0.000 

0.000 

0.2 

0.339 

0.000 

0.000 

0.3 

0.327 

0.000 

0.000 

0.4 

0.325 

0.000 

0.000 

0.5 

0.331 

0.000 

0.000 

W/NiFe  Type-II 

0.1 

0.340 

1.125 

5.000 

0.2 

0.314 

1.250 

5.000 

0.3 

0.304 

1.375 

5.000 

0.4 

0.298 

4.750 

2.875 

0.5 

0.292 

4.875 

2.875 

W/4340  Steel  Type-I 

0.3 

0.364 

1.125 

5.000 

W/Gu  Type-I 

0.3 

0.218 

0.625 

5.000 

W/Fe  Type-I 

0.3 

0.268 

1.000 

5.000 

W,  r  <  5 

NiFe,  r  >  5 

— 

0.194 

2.875 

4.000 

W,  r  >  5 

NiFe,  r  <  5 

— 

0.154 

4.875 

1.500 

Of  the  FGMs  studied,  the  type-I  with  NiFe  matrix  delays  most  the  onset  of  an  ASB  and 
that  with  Fe  matrix  the  least.  Figure  3.6  depicts  contours  of  the  effective  plastic  strain  in 
the  type-I  FGMs  with  the  matrix  material  either  4340  steel,  OFHG  copper,  or  Fe.  In  each 
case  an  ASB  initiates  from  a  different  point  on  the  top  surface. 

Gomposites  comprised  of  W  in  the  region  r  <  5  mm  and  NiFe  elsewhere  or  vice-a-versa 
were  also  analyzed.  ASBs  initiated  in  these  two  composites  at  much  lower  values  of  the 
average  axial  strain  than  in  all  of  the  other  FGMs  and  in  the  two  homogeneous  bodies. 
Recall  that  material  properties  are  discontinuous  across  the  surface  r  =  H.  It  suggests  that 
the  presence  of  a  material  singular  surface  enhances  the  initiation  of  an  ASB.  The  gradients 
of  displacements  and  temperature  must  be  discontinuous  across  r  =  H  in  order  for  surface 
tractions  and  the  normal  component  of  the  heat  flux  to  be  continuous. 

The  significant  difference  between  the  initiation  times  of  ASBs  in  the  two  composites  with 
W  in  the  region  r  <  5  mm  and  NiFe  in  r  >  5  mm  or  vice-versa  can  not  be  explained  on  the 
basis  of  the  arrival  time  of  the  elastic  loading  wave  since  \JE/ p  is  nearly  the  same  for  W  and 
NiFe.  However,  their  acoustic  impedances,  \/E~p,  are  quite  different.  Batra  and  Kwon  [88] 
scrutinized  the  initiation  and  development  of  ASBs  in  simple  shearing  deformations  of  a 
bimetallic  body  with  the  lower  half  made  of  one  thermoviscoplastic  material  and  the  upper 
half  of  a  different  thermoviscoplastic  material.  For  both  materials  having  the  same  value  of 
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the  shear  modulus  and  the  mass  density,  they  found  that  the  ASB  developed  in  the  body 
with  the  higher  value  of  the  thermal  softening  coefficient  even  though  the  initial  nonuniform 
temperature  was  symmetric  about  the  centerline  and  maximum  at  the  center.  The  differences 
in  thermal  conductivities  shifted  the  ASB  a  little.  No  shift  in  the  center  of  the  ASB  away 
from  the  centerline  occurred  even  when  the  shear  modulus  of  one  material  equaled  twice 
that  of  the  other  material.  The  data  included  in  the  last  two  rows  of  Table  3.4  indicates 
that  an  ASB  initiated  at  the  point  (2.875,  4.0)  or  r  =  4.93  mm  when  W  occupied  the  region 
r  <  5  mm,  and  at  the  point  (4.875,  1.50)  or  r  =  5.10  mm  for  W  in  the  region  r  >  5  mm. 
The  coordinates  of  points  are  in  the  reference  conhguration.  These  results  agree  with  those 
of  Batra  and  Kwon  [88].  Zhou  et  al.  [89]  analyzed  numerically  and  experimentally  plane 
strain  deformations  of  a  particulate  composite  comprised  of  W  particles  in  NiFe  matrix  and 
modeled  each  constituent  as  a  thermoviscoplastic  material.  Thus  there  were  several  surfaces 
of  discontinuity  present  in  the  body.  They  also  found  that  an  ASB  initiated  sooner  in  the 
composite  than  in  either  one  of  the  two  constituents.  They  neither  homogenized  material 
properties  nor  investigated  the  effect  of  the  spatial  variation  of  the  volume  fraction  of  the 
two  constituents  nor  considered  the  evolution  of  porosity. 

An  ASB  initiates  from  the  point  at  the  center  of  the  shaded  area  in  Figure  3.7a.  The 
shaded  zones  in  Figure  3.7b-d  depict  the  shear  banded  regions.  Since  the  length  of  this 
region  increases  in  both  directions  around  the  point  shown  in  Figure  3.7a,  one  can  conclude 
that  the  shear  band  propagates  outwards  in  both  directions. 

For  each  problem  studied,  the  time  step  was  seen  to  drop  drastically  once  an  ASB  had 
initiated  at  a  point  in  the  body.  The  numerical  values  of  the  ASB  initiation  time  may  depend 
upon  the  FE  mesh  employed.  However,  differences  in  results  from  a  much  hner  mesh  than 
the  one  used  here  are  likely  to  be  less  than  5%.  Recall  that  for  a  type-I  FGM,  as  the  ASB 
initiation  time  obtained  with  6400  uniform  elements  differed  from  that  with  1600  uniform 
elements  by  less  than  2.2%.  One  could  use  either  adaptively  rehned  meshes,  e.g.  see  [21],  or 
a  meshless  method  such  as  the  modihed  smoothed  particle  hydrodynamics  [90]  to  compute 
results.  The  former  technique  smoothens  deformation  helds  and  tends  to  delay  the  ASB 
initiation  but  gives  narrower  bands.  The  latter  is  computationally  less  efficient  than  the 
FEM. 

Other  compositional  prohles  symmetric  about  the  two  centroidal  axes  include  Vf^NiPe  = 
a  +  where  a,  b,  p  and  q  are  real  numbers  such  that  0  <  Pf^NiPe  <  1  for  all  admissible 

values  of  Xi  and  X2.  These  have  not  been  examined. 
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Figure  3.4:  Contours  of  the  effective  plastic  strain  in  type-I  FGMs  with  the  maximum 
volume  fraction  of  NiFe  at  r  =  H  equal  to  (a)  0.1,  (b)  0.2,  (c)  0.3  (d)  0.4, 
and  (e)  0.5 
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Figure  3.6:  Contours  of  the  effective  plastic  strain  in  type-I  FGMs  with  the  maximum 
volume  fraction  of  W  at  r  =  0  equal  to  0.7  and  (a)  4340  steel,  (b)  Cu,  and 
(c)  Fe  as  the  other  constituent 
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(a) 


(b) 


Figure  3.7:  Shear-banded  region,  shown  shaded,  in  a  W/NiFe  composite  with  W  in  the 
region  r  <  5  mm  and  NiFe  in  r  >  5  mm  at  (a)  42/^s,  (b)  46/^s,  (c)  50^s,  (d) 
54^s,  and  (e)  58/is;  the  effective  plastic  strain  in  the  shaded  region  exceeds 
1.0 
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3.4  Conclusions 


We  have  analyzed  the  initiation  and  propagation  of  ASBs  in  FG  bodies  deformed  in  plane 
strain  tension  at  an  average  axial  strain-rate  of  5000/s.  Each  FGM  studied  is  comprised  of 
two  constituents  one  of  which  is  tungsten  and  the  other  NiFe,  copper,  4340  steel  or  iron. 
Several  spatial  distributions  of  the  two  constituents  in  the  region  r  <  H  are  considered. 
Here  r  is  the  radial  distance  from  the  centroid  of  the  square  cross-section  of  side  2H.  The 
material  is  homogeneous  in  the  region  r  >  H.  For  uniform  40  x  40  and  120  x  120  hnite 
element  meshes  of  4-node  quadrilateral  elements,  the  ASB  initiation  time  differed  by  2.1%. 
It  is  found  that  an  ASB  initiates  either  at  the  specimen  centroid  or  at  a  point  on  the  top 
surface  where  the  normal  component  of  velocity  is  prescribed.  In  each  case,  it  propagates 
along  the  direction  of  the  maximum  shear  stress.  An  ASB  initiates  in  a  FGM  at  a  lower 
value  of  the  average  axial  strain  than  that  in  a  homogeneous  body  made  of  either  one  of  the 
two  constituents.  Thus,  for  the  material  combinations  studied  here,  the  addition  to  tungsten 
of  a  material  that  is  less  susceptible  to  adiabatic  shear  banding  than  tungsten,  enhances  the 
susceptibility  of  the  FGM  to  adiabatic  shear  banding.  The  analysis  can  thus  be  used  to 
optimize  the  compositional  prohle  to  either  delay  or  to  enhance  the  onset  of  an  ASB. 
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Chapter  4 

Brittle  and  Ductile  Fracture  in 

•  *1 

Functionally  Graded  Materials 

4.1  Problem  Statement 

Here  we  investigate  the  initiation  and  propagation  of  brittle  and  ductile  failures  in  tungsten/ 
nickel-iron  FGMs.  The  failure  of  FGMs  is  of  great  interest  to  designers-FGMs  allow  material 
transitions  without  a  failure-prone  interface,  but  little  is  known  about  crack  propagation  and 
material  failure  of  the  FGMs  themselves.  Both  constituents  and  the  overall  composite  are 
taken  to  be  thermoelastoviscoplastic,  allowing  the  plasticity-driven  ductile  crack  to  initiate 
and  propagate.  Very  few  researchers  have  simulated  it  numerically.  The  crack  propagation 
speed  is  obtained  in  both  pure  materials  and  in  the  FGM,  and  the  interaction  between  ASBs 
and  the  ensuing  ductile  fracture  is  examined. 

4.2  Numerical  Tests 

4.2.1  Initial  and  Boundary  Conditions 

The  plane  strain  FE  code  described  in  Section  1.2  is  utilized;  the  crack  propagation  routines 
and  the  contact  algorithm  described  in  Section  1.2.6  are  implemented.  For  both  brittle  and 
ductile  fracture  tests,  the  body  is  assumed  to  be  initially  at  rest,  stress-free,  at  a  uniform 
temperature,  and  have  no  initial  porosity.  Under  the  plane  strain  assumption,  it  is  tacitly 
assumed  that  the  volume  fraction  of  constituents  is  constant  in  the  X3  direction.  For  the 
brittle  fracture  case,  we  examine  a  center-cracked  square  cross-section  of  length  and  width 

^Material  in  this  Chapter  is  from  Ref.  [91]. 
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Figure  4.1:  (a)  FE  mesh  for  plane  strain  tensile  deformations  of  a  pre-cracked  plate;  inset 
shows  details  of  the  mesh  around  the  cracked  region;  (b)  FE  mesh  for  plane 
strain  shear  deformations  of  a  plate 


2H  deformed  in  plane-strain  tension.  Symmetry  along  the  X2  axis  is  exploited;  deformations 
in  half  of  the  body,  shown  in  Figure  4.1a,  are  analyzed.  The  symmetry  about  the  Xi  axis  is 
not  exploited  due  to  the  nature  of  the  crack  propagation  algorithm,  which  only  allows  cracks 
internal  to  the  body.  Boundary  conditions  arising  from  the  symmetry  of  deformations  are 
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imposed  at  points  on  the  centroidal  axis  Xi  =  0.  The  other  vertical  surface  Xi  =  H  is  taken 
to  be  traction  free  and  thermally  insulated.  Normal  velocity,  null  tangential  tractions  and 
zero  heat  flux  are  prescribed  on  the  top  horizontal  surface  X2  =  H.  The  prescribed  normal 
velocity  increases  linearly  with  time  to  its  steady  state  value  Vq  in  1  fis  and  is  then  held 
hxed.  Thus,  for  the  brittle  fracture  problem 


T21  = 
T21  = 

Tu  = 


Tn  = 

0,  vi 

0,  Q2 


0,  Q 

=  0, 

=  0, 


1  =  0  on  Xi  =  iJ, 

Qi  =  0  on  Xi  =  0, 

f  vot,  0  <t  <  1  /is, 
V2=  < 

Vo,  t>  IflS, 


on  X2 


H. 


(4.1) 


For  the  ductile  fracture  case,  we  analyze  plane  strain  shear  deformations  of  a  body  of  length 
2L  and  height  2H,  as  shown  in  Figure  4.1b.  Periodic  boundary  conditions  are  applied  on 
the  surfaces  Xi  =  ±L,  which  specihes  that  the  displacements  and  temperatures  are  equal  at 
corresponding  points-this  simulates  an  inhnitely  long  specimen  of  height  2H.  The  top  and 
bottom  surfaces  are  taken  to  be  thermally  insulated  and  restrained  from  motion  in  the  X2 
direction.  Equal  and  opposite  tangential  velocity  vi  that  increases  linearly  with  time  from 
zero  to  its  steady  value  Vq  is  1/is  is  applied.  Thus  on  X2  =  X2  =  ±hf,  we  have 


Q2  =  0,  V2  =  0,  Vi 


±Vot,  0  <  t  <  1  IIS, 
Uo,  t  >  1  IIS. 


(4.2) 


For  t  >  1  IIS,  the  average  shear  strain-rate  is  vq/H.  Note  that  only  the  top  and  the  bottom 
surfaces  are  restrained  from  moving  vertically  and  other  material  particles  are  free  to  move 
in  the  a;ia;2-plane. 

As  in  Chapter  3,  the  local  material  properties  were  derived  by  the  rule  of  mixtures. 
As  shown  in  Chapter  2,  this  provides  a  reasonable  estimate  of  the  thermoelastoviscoplastic 
properties  over  a  wide  range  of  volume  fractions,  and  the  properties  are  easily  computed 
from  the  properties  of  the  constituents. 


4.2.2  Failure  Initiation  Criterion 

Brittle  failure 

Ritchie  et  al.  [92]  proposed  that  the  brittle  failure  initiates  at  a  point  when  ap/ao  =  3.0  over 
a  certain  length  that  depends  upon  the  microstructure  of  the  material  and  usually  equals 
a  grain  diameter.  Here  ap  and  ao  equal,  respectively,  the  maximum  principal  tensile  stress 
at  a  point  and  the  yield  stress  in  a  quasistatic  simple  tension/compression  test.  For  the 
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Johnson-Cook  material,  (Tq  can  be  taken  to  equal  A.  Tensile  experiments  of  Hendrickson  et 
al.  [93]  on  prenotched  steel  plates  with  a  yield  stress  of  705  MPa  deformed  at  nominal  stress 
rates  of  ~  1  to  10^  MPa/s  revealed  that  brittle  failure  occurred  when  (Jp/ao  —  2.34.  This 
value  of  (Tp/ (To  was  found  to  be  independent  of  the  temperature  and  the  rate  of  loading.  Here 
we  assume  that  brittle  failure  initiates  at  a  point  when  ap/a^  =  3.0  there  and  propagates  in 
the  direction  of  the  minimum  gradient  in  (jp. 

Ductile  failure 

Ductile  failure  is  generally  believed  to  initiate  due  to  the  nucleation  and  coalescence  of  voids 
within  an  ASB.  We  recall  that  the  effective  plastic  strain  induced  within  an  ASB  exceeds  1. 
Accordingly,  a  ductile  failure  is  assumed  to  ensue  at  a  point  when  the  effective  plastic  strain 
equals  1.5  and  propagate  in  the  direction  of  the  minimum  gradient  in  the  effective  plastic 
strain. 


4.3  Computation  and  Discussion  of  Results 

4.3.1  Wave  Propagation  in  a  Linear  Elastic  EG  Bar 

We  validate  our  methodology  of  analyzing  transient  deformations  of  a  FG  body  by  studying 
wave  propagation  in  a  linear  elastic  bar  whose  Young’s  modulus,  E,  mass  density  p,  and 
Poisson’s  ratio  u  are  given  by 


E  =  200(1  +  0.25Xi)GPa,  p  =  10^(1  +  0.25Xi)-ikg/m^  z/  =  0.29. 


(4.3) 


In  equation  (4.3)  Xi  is  in  mm.  Thus  the  wave  speed  Cw  = 


E{l-iy) 

; - - - -  varies 

y(1  +  '^)(1  -  2z^)/ 

affinely  from  5.12  km/s  at  Xi  =  0  to  6.40  km/s  at  Xi  =  20  mm.  However,  the  acoustic 
impedance  (=  \/Ep)  is  constant  throughout  the  bar.  In  order  to  simulate  one-dimensional 
elastic  deformations  of  the  bar,  all  nodes  were  restrained  from  moving  in  the  X2-direction, 
and  the  material  parameter  A  was  assigned  a  very  high  value  to  suppress  plastic  deformations. 
Since  the  computer  code  for  plane  strain  deformations  tacitly  assumes  that  U3  =  0,  therefore, 
the  only  non-vanishing  component  of  displacement  is  Ui.  The  20  mm  x  0.5  mm  bar  was 
divided  into  800  x  20  uniform  4-node  quadrilateral  elements.  Various  integrals  appearing  in 
the  weak  formulation  of  the  problem  were  computed  by  using  the  2x2  quadrature  rule  with 
values  of  material  parameters  evaluated  at  the  integration  points.  The  same  strategy  was 
employed  by  Batra  [94]  to  analyze  static  hnite  plane  strain  deformations  of  an  inhomogeneous 
Mooney-Rivlin  material.  An  analytical  solution  of  the  problem,  using  the  Laplace  transform 
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axial  velocity  (km/s) 


Figure  4.2:  (a)  Time  history  of  the  axial  velocity  prescribed  at  the  left  end  of  the  bar;  (b) 
comparison  of  the  time  histories  of  the  axial  stress  at  x  =  10  mm  obtained 
from  the  numerical  and  analytical  solutions;  (c)  comparison  of  the  spatial 
variation  in  the  wave  speed 
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technique,  has  been  given  by  Chiu  and  Erdogan  [95].  It  is  inverted  numerically  by  employing 
Laguerre  polynomials  to  hnd  the  axial  stress  and  the  axial  velocity  as  a  function  of  time  t. 

Figure  4.2a  depicts  the  axial  velocity  prescribed  at  the  end  Xi  =  0  of  the  bar.  The  time 
histories  of  the  axial  stress  at  Xi  =  10  mm  obtained  from  the  numerical  and  the  analytical 
solutions  are  compared  in  Figure  4.2b,  and  the  spatial  variations  of  the  wave  speed  obtained 
from  the  numerical  and  the  analytical  solutions  are  compared  in  Figure  4.2c.  Except  for  small 
oscillations,  the  two  sets  of  results  agree  well  with  each  other.  The  amplitude  of  oscillations 
can  be  diminished  by  introducing  an  artihcial  viscosity,  but  this  was  not  attempted. 

4.3.2  Simulation  of  Brittle  Fracture 

Crack  Propagation  Speed 

The  10  mm  x  10  mm  cross-section  has  an  initial  sharp  crack  of  length  1  mm  at  the  horizontal 
centroidal  axis  with  the  center  of  the  crack  coincident  with  the  centroid  of  the  plate.  The 
plate  is  deformed  at  an  average  axial  strain-rate  of  either  200/s  or  2,000/s.  Values  of  Young’s 
modulus,  E,  mass  density,  p,  and  Poisson’s  ratio,  z/,  for  W  and  NiFe,  and  the  speeds  of  the 
longitudinal  wave  in  a  bar  and  the  Rayleigh  wave  speed  are  listed  in  Table  4.1;  values  of 
other  material  parameters  are  given  in  Table  4.2  and  in  equation  (4.4).  We  note  that  NiFe 
exhibits  considerably  higher  strain-  and  strain-rate  hardening  than  W. 


Table  4.1:  Material  parameters  and  wave  speeds  for  nickel-iron  and  tungsten 


Material 

Young’s 

modulus 

(GPa) 

Poisson’s 

ratio 

Mass 

density 

(kg/m^) 

Bar  wave 
speed 
(m/s) 

Acoustic 
impedance 
(kg/m^  s) 

Rayleigh 
wave  speed 
(m/s) 

NiFe 

255 

0.29 

9200 

5265 

48.44  X  10^ 

3035 

Tungsten 

400 

0.29 

19300 

4552 

87.86  X  10^ 

2624 

Table  4.2:  Values  of  material  parameters  for  NiFe  and  W 


Material 

A 

(MPa) 

B 

(MPa) 

C 

n 

(K) 

c 

(J/kg  K) 

K 

(W/mK) 

a 

(lO'VK) 

m 

n 

NiFe 

150.0 

546.0 

0.0838 

1225 

382 

100 

15 

1.0 

0.208 

W 

730 

562 

0.290 

1723 

138 

160 

5.3 

1.0 

0.075 
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Other  parameters  were  assigned  the  following  values  for  both  materials. 

eg  =  IQ-Vs,  (3i  =  1.5,  (52  =  1.0,  /2  =  0.04,  S2  =  0.1,  Or  =  273  K, 

(4.4) 

r  =  lO-^s,  =  0.5,  /,  =  0.15,  U  =  2/3,  fj  =  0.25,  m  =  1.0. 

The  FE  mesh  used  to  analyze  the  problem  is  depicted  in  Figure  4.1a,  details  of  mesh 
around  the  crack-tip  are  given  in  the  inset.  The  mesh  consists  of  17,444  4- node  quadrilateral 
elements  with  1,080  elements  along  the  axis  of  the  crack.  There  are  108  uniform  elements 
behind  the  crack-tip  and  972  ahead  of  it;  thus  the  length  of  an  element  is  4.63  /rm.  The 
appropriateness  of  the  mesh  has  been  ascertained  by  ensuring  that  the  computed  speed  of 
an  elastic  wave  is  very  close  to  the  analytic  value  when  E  and  p  vary  along  the  direction  of 
propagation  of  the  wave.  Once  the  brittle  failure  criterion  at  a  node  is  met,  that  node  is  split 
into  two  essentially  overlapping  but  unconnected  nodes  as  described  in  Section  1.2.6.  Thus  an 
elastic  unloading  wave  emanates  from  the  newly  created  crack  surfaces  and  propagates  into 
the  body.  The  position  of  the  crack-tip  at  different  times  is  determined  and  a  polynomial  is  £t 
to  the  data.  The  hrst  derivative  of  this  polynomial  fit  gives  speed,  of  crack  propagation  as 
a  function  of  time  or  the  position  of  the  crack-tip.  It  was  found  that  Cs  so  determined  is  very 
sensitive  to  the  polynomial  fit.  Three  curve  fits  with  the  coefficient  of  regression  >  0.999  gave 
noticeably  different  values  of  Cg-  Thus  Cg  at  a  point  is  taken  to  equal  the  slope  of  the  least 
squares  line  through  the  neighboring  ±10  points.  Plates  consisting  of  pure  W  and  pure  NiFe 
were  tested  along  with  FG  plates  consisting  of  100%  W  at  Xi  =  0  to  0%  W  at  Xi  =  5mm 
(hereafter  referred  to  as  W2NiFe)  and  vice-versa  (hereafter  referred  to  as  NiFe2W).  The 
Rayleigh  wave  speed  and  the  crack  propagation  speed,  Cg,  versus  the  crack  length  for  nominal 
strain-rates  of  200/s  and  2,000/s  are  depicted  in  Figure  4.3a-d.  Freund’s  [96]  analysis  of  crack 
propagation  in  an  inhnite  homogeneous  linear  elastic  body  shows  that  the  maximum  crack 
propagation  speed  equals  the  Rayleigh  wave  speed;  Eischen  [49]  has  proved  a  similar  result 
for  inhomogeneous  materials.  For  each  of  the  homogeneous  and  the  two  FG  plates  Cg  for 
the  nominal  strain-rate  of  2,000/s  is  higher  than  that  for  the  nominal  strain-rate  of  200/s. 
For  the  W  plate  deformed  at  200/s,  Cg  increases  as  the  crack  propagates  to  the  right  edge 
of  the  plate  but  at  the  higher  strain-rate  of  2,000/s,  it  soon  approaches  a  steady  value  that 
is  a  little  less  than  the  Rayleigh  wave  speed.  For  a  NiFe  plate  deformed  at  2,000/s,  the 
crack  begins  to  propagate  to  the  right  and  then  a  large  region  of  the  plate  ahead  of  the 
crack  fails  instantaneously  as  indicated  by  the  maximum  principal  tensile  stress  exceeding 
3(Jo  simultaneously  everywhere  in  this  region.  It  is  signified  in  Figure  4.3b  by  the  sudden 
drop  in  the  crack  propagation  speed  Cg]  the  shattered  region  is  depicted  in  Figure  4.3e  as 
unshaded. 
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crack  propagation  speed  (km/s)  crack  propagation  speed  (km/s) 


Figure  4.3:  (a)-(d):  Crack  propagation 
speed  versus  crack  length  at  nominal  strain- 
rates  of  200/s  and  2000/s  for  (a)  pure  tung¬ 
sten,  (b)  pure  nickel-iron,  (c)  W2NiFe,  and 

(d)  NiFe2W  (key: - Rayleigh  wave  speed; 

- -  nominal  strain-rate  =  2000/s;  -  •  - 

•  -  nominal  strain-rate  =  200/s);  (e)  shat¬ 
tered  region  in  the  NiFe  plate  deformed  at 
a  nominal  strain-rate  of  2000 /s  is  unshaded 


(e) 
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For  the  W2NiFe  FG  plate,  Cg  continues  to  increase  with  the  crack  extension,  is  always 
less  than  the  Rayleigh  wave  speed  when  the  plate  is  deformed  at  200/s  but  approaches 
the  Rayleigh  wave  speed  when  the  nominal  strain-rate  is  2,000/s.  Except  for  differences  in 
magnitudes,  the  curves  in  Figures  4.3a  and  4.3c  are  similar.  Thus  as  far  as  crack  propagation 
due  to  brittle  failure  is  concerned,  the  W  and  the  W2NiFe  FG  plates  behave  alike.  However, 
for  the  NiFe2W  FG  plate,  even  though  the  Rayleigh  wave  speed  decreases  monotonically  with 
the  distance  from  the  left  edge  because  of  the  spatial  variation  in  the  material  properties, 
the  computed  crack  speed  Cg  first  increases  and  approaches  essentially  a  steady  value  after 
the  crack  has  propagated  a  certain  distance.  The  crack  speed  is  higher  when  the  nominal 
strain-rate  equals  2,000/s  than  that  when  the  nominal  strain-rate  is  200/s.  No  shattering  of 
the  NiFe2W  FG  plate  at  either  of  the  two  strain-rates  is  observed. 

Shattering  Phenomenon 

Returning  to  the  shattering  of  the  NiFe  plate,  several  numerical  experiments  were  conducted 
by  varying  the  nominal  axial  strain-rate.  The  plate  shattered  at  nominal  axial  strain-rates 
exceeding  1,130/s  but  did  not  shatter  at  a  nominal  axial  strain-rate  of  1,120/s.  The  crack 
extension/elongation  at  the  instant  of  the  plate  shattering  decreased  with  an  increase  in  the 
nominal  axial  strain-rate;  it  equaled  ~  2.01,  1.47  and  0.741  mm  for  nominal  axial  strain-rates 
of  1,130,  2,000  and  3,000/s,  respectively. 

Axial  Force 

We  have  plotted  in  Figure  4.4  the  variation  of  the  computed  axial  load  as  the  crack  propagates 
to  the  right.  Small  oscillations  in  the  load  due  to  the  arrival  of  unloading  waves  emanating 
from  the  newly  formed  crack  surfaces  have  been  diminished  by  using  a  neighboring-point- 
average  technique,  included  in  TecPlot.  At  a  nominal  strain-rate  of  200/s,  the  crack  prop¬ 
agation  is  stable  in  both  pure  W  and  W2NiFe  FG  plates  as  signihed  by  either  the  load 
remaining  essentially  steady  or  increasing  except  when  the  crack  has  propagated  to  a  point 
near  the  right  edge.  At  the  higher  strain-rate  of  2,000/s,  the  axial  load  for  pure  W  remains 
nearly  constant  at  ~  4.67  kN  until  the  crack-tip  has  advanced  by  1.75  mm  and  then  slowly 
increases  until  the  crack-tip  is  close  to  the  outer  edge  of  the  plate.  When  the  precracked  W 
plate  is  pulled  at  a  nominal  axial  strain-rate  of  200/s,  the  axial  load  continues  to  increase 
slowly  until  the  crack  length  equals  2.2  mm  and  then  gradually  decreases.  Recall  that  the 
surface  area  of  the  horizontal  plane  containing  the  crack  decreases  as  the  crack  extends,  and 
except  for  the  hrst  1  /xs,  the  bar  is  being  pulled  with  a  uniform  axial  velocity.  Thus  the 
working  of  external  forces  is  proportional  to  the  axial  force.  The  variation  of  the  axial  load 
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with  the  advance  of  the  crack  in  NiFe  and  NiFe2W  FG  plates  is  similar  to  that  in  the  W 
plate  except  that  a  large  portion  of  material  ahead  of  the  crack  suddenly  fails  in  the  pure 
NiFe  plate  deformed  at  a  nominal  strain-rate  of  2,000/s.  For  the  NiFe2W  FG  plate,  the 
variation  of  the  axial  load  with  the  crack  length  is  similar  to  that  for  the  W2NiFe  FG  plate. 


(b) 


(d) 


Figure  4.4:  Axial  load  versus  crack  length  at  nominal  strain-rates  of  200/s  and  2000/s 
for  (a)  tungsten,  (b)  nickel- iron,  (c)  W2NiFe,  and  (d)  NiFe2W  plates 


J-integral 

A  detailed  examination  of  stresses  and  strains  induced  within  the  specimen  revealed  that 
deformations  were  virtually  elastic  everywhere  except  near  the  crack-tip  where  small  plastic 
strains  developed.  In  order  to  see  if  concepts  of  the  strain  energy  release  rate  or  the  J-integral 
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(b) 


Figure  4.5:  Variation  of  the  J-integral  with  the  crack  extension  in  a  W  plate 


can  be  used,  we  computed  the  J-integral  defined  by 


J  = 

w  = 


lim  /  {w  +  T)dx2  —  (TijUj^^^ds^ 


r^o 

pt 


dxi 


(4.5) 


Here  F  is  a  closed  contour  enclosing  the  crack-tip,  and  ds  is  an  element  of  arc  length  on 
r.  Figure  4.5a  shows  the  value  of  the  J-integral  as  a  function  of  the  crack-tip  location  for 
the  pure  W  specimen  deformed  at  nominal  strain-rates  of  200  and  2,000/s.  The  contour  F 
used  to  evaluate  the  J-integral  appeared  stationary  to  an  observer  situated  at  the  crack- 
tip  and,  hence,  moving  with  it.  The  contour  F  spanned  5  elements  behind  the  crack-tip, 
5  elements  ahead  of  it  and  10  elements  perpendicular  to  the  crack-axis  (each  element  is  a 
square  of  side  length  4.63  /rm  and  the  starter  crack  of  length  1000  /rm).  For  an  average  strain- 
rate  of  200/s,  the  J-integral  increases  monotonically  signifying  an  increasing  resistance  to 
crack  propagation.  However,  at  the  average  axial  strain-rate  of  2,000/s,  the  J-integral  first 
increases  with  an  advance  in  the  crack-tip  but  then  oscillates  wildly.  These  oscillations  can 
not  be  attributed  to  the  contour  F  used  to  evaluate  J  since  a  similar  contour  with  respect 
to  the  crack-tip  gave  reasonable  values  of  the  J  integral  for  all  locations  of  the  crack-tip  at 
the  lower  strain-rate  of  200/s  and  for  earlier  stages  of  the  crack  propagation  in  the  plate 
deformed  at  an  axial  strain-rate  of  2,000/s.  The  time  history  of  the  evolution  of  the  J- 
integral  is  exhibited  in  Figure  4.5b  where  its  values  at  the  instants  of  crack  initiation  at 
different  points  are  shown.  Values  of  the  J-integral  at  the  instant  of  the  first  crack  initiation 
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in  the  W-plate  are  the  same  and  equal  ~  0.22  kJ/m^  at  nominal  strain-rates  of  200  and 
2,000/s. 

4.3.3  Simulation  of  Ductile  Failure 

In  order  to  simulate  crack  propagation  due  to  ductile  failure  we  analyze  plane  strain  thermo¬ 
mechanical  deformations  of  a  block  of  material  with  equal  and  opposite  tangential  velocities 
prescribed  on  its  top  and  bottom  surfaces  as  described  in  Section  4.2.1.  Thus  the  steady- 
state  nominal  strain-rate  is  Vq/H,  and  equals  5,000/s  for  simulations  discussed  herein.  The 
FE  mesh  used  to  analyze  the  problem  is  depicted  in  Figure  4.1b.  The  elements  are  of  size 
4.6  X  4.6  /xm  in  the  central  and  end  portions  but  get  larger  in  other  regions.  For  the  FG 
bodies,  the  composition  varies  in  the  horizontal  direction  only  from  either  0%  W  at  the 
centroid  to  100%  W  at  the  edges  or  vice-versa.  The  yield  stress  of  the  material  in  an  area 
4.6  /xm  thick  and  200  /rm  long  located  symmetrically  about  the  two  centroidal  axes  was 
reduced  by  30%;  this  defect  acts  as  a  nucleation  spot  for  an  ASB.  Recall  that  a  crack  is 
assumed  to  open  at  a  point  when  the  effective  plastic  strain  there  equals  1.5;  thus  the  crack 
opens  inside  the  body,  generating  two  crack-tips.  It  was  found  that  the  ASB-tips  and  the 
crack-tips  propagate  horizontally  with  equal  and  opposite  velocities;  thus  their  propagation 
to  the  right  is  described  below. 

Crack  speed 

Figures  4.6  and  4.7  depict  the  variation  in  the  crack  length,  the  ASB  length,  the  crack-tip 
speed  and  the  ASB  speed  as  they  advance  in  the  two  homogeneous  materials,  namely  W 
and  NiFe,  and  the  FGM  W2NiFe.  For  the  NiFe2W  FG  plate,  in  spite  of  the  rather  strong 
defect  at  the  centroid,  the  crack  originated  from  points  (— T,  0)  and  (L,  0)  and  propagated 
inwards.  Thus  it  propagated  from  W  rich  region  to  NiFe  rich  region,  and  the  situation 
is  similar  to  that  for  crack  propagation  in  the  W2NiFe  FG  plate.  An  ASB  initiates  at  a 
point  when  the  effective  plastic  strain  there  equals  1.5;  thus  the  difference  between  an  ASB 
and  a  crack  is  that  new  thermally  insulated  traction  free  surfaces  are  created  when  a  crack 
opens  but  the  material  is  intact  during  the  initiation  and  propagation  of  an  ASB.  The  crack 
surfaces  are  always  in  contact  with  one  another  due  to  a  compressive  normal  stress  acting 
on  them.  In  W,  the  time  histories  of  the  crack  length  and  the  ASB  length  and  hence  their 
speeds  of  propagation  are  virtually  identical.  However,  in  NiFe  the  crack  length  is  larger 
than  the  ASB  length.  Hence  the  crack  propagation  speed  is  higher  than  the  ASB  speed. 
Except  when  the  ASB/crack  arrives  at  xi  =  ±L,  its  speed  increases  from  ~  0.1  km/s  at 
origination  to  ~  1.8  km/s  at  the  end.  The  crack  accelerates  slowly  in  the  beginning  but 


60 


(a) 


(b) 


Figure  4.6:  Time  histories  of  the  evo¬ 
lution  of  the  ASB  and  the  crack  in  (a) 
W,  (b)  NiFe  and  (c)  W2NiFe  plates  de¬ 
formed  in  plane  strain  shear 


quite  rapidly  towards  the  end.  Note  that  an  ASB/crack  initiates  in  NiFe  at  165  /is  but  in  W 
at  35  JUS.  Interestingly,  the  time  histories  of  the  crack  propagation  speeds  in  W  and  NiFe  are 
nearly  coincident.  Whether  it  is  so  in  all  homogeneous  materials  remains  to  be  seen.  Had 
we  assumed  that  an  ASB  forms  at  the  effective  plastic  strain  of  1.0  and  a  crack  opens  at  the 
effective  plastic  strain  of  1.5,  the  ASB  length  would  have  been  larger  than  the  crack  length. 
However,  their  speeds  of  propagation  would  be  essentially  the  same. 

In  W2NiFe  FG  plate,  the  ASB/crack  initiation  time  is  between  those  for  the  pure  W 
and  the  pure  NiFe  plates;  a  similar  situation  occurs  in  the  brittle  failnre  of  these  materials. 
The  ASB/crack  propagation  speed  in  the  FG  plate  is  signihcantly  less  (approximately  by  an 
order  of  magnitude)  than  that  in  either  of  its  two  constituents. 
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We  note  that  both  the  crack-tip  and  the  ASB-tip  do  not  advance  monotonically.  Rather, 
once  a  crack  has  opened  at  a  node,  it  takes  a  little  while  for  the  fracture  criterion  to  be  met 
at  the  node  immediately  ahead  of  the  crack-tip.  Thus,  the  crack  length  as  a  function  of  time 
exhibits  a  stair-case  pattern.  Results  plotted  in  Figure  4.6  have  been  smoothened. 

For  ductile  fracture,  there  are  signihcant  plastic  deformations  induced,  so  the  J-integral 
is  not  evaluated. 


Tangential  force 

We  have  plotted  in  Figure  4.8  time-histories  of  the  tangential  force  required  to  deform  the 
specimen.  Since  the  top  and  the  bottom  surfaces  are  constrained  from  moving  in  the  vertical 
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(b) 


Figure  4.8:  Time  histories  of  the  tan¬ 
gential  force  in  (a)  W,  (b)  NiFe  and  (c) 
W2NiFe  plates  deformed  in  plane  strain 
shear 


(c) 


direction,  the  working  of  external  forces  is  due  to  tangential  tractions.  Furthermore,  this 
working  is  proportional  to  the  tangential  force  as  the  tangential  velocity,  except  for  the  hrst  1 
/is,  is  constant.  As  expected  the  tangential  force  decreases  with  the  opening  of  a  crack  and  it 
continues  to  decrease  as  the  crack  elongates.  With  the  crack  extension,  smaller  surface  area 
supports  the  external  load.  Also,  with  continued  plastic  deformations,  the  material  softens 
and  its  capacity  to  support  external  load  diminishes.  The  rate  of  decrease  of  the  tangential 
force  is  higher  for  W  than  that  for  NiFe  since  the  strain  and  strain-rate  hardening  effects  are 
higher  in  NiFe  than  in  W.  With  the  extension  of  the  crack,  the  tangential  force  drops  more 
rapidly  for  the  FG  plate  than  that  for  the  two  homogeneous  plates. 
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Figure  4.9:  Effective  strain-rates  ahead  of  the  tip  in  plates  (a)  with  crack  and  (b)  with 
ASB  only;  here  a  is  the  current  length  of  the  crack 


4.3.4  Effect  of  Crack  Opening  on  Deformation  Fields 

In  order  to  delineate  the  differences,  if  any,  between  the  deformation  fields  ahead  of  the  ASB- 
tip  and  the  crack-tip,  we  have  plotted  in  Figure  4.9  the  spatial  variation  of  the  effective  plastic 
strain-rate  in  a  W  plate  deformed  at  a  nominal  strain-rate  of  5000/s  with  and  without  crack 
propagation.  Because  of  the  unloading  elastic  waves  emanating  from  newly  formed  crack 
surfaces,  the  strain-rates  ahead  of  the  crack-tip  are  oscillatory;  cf.  Figure  4.9a.  The  slope 
of  the  line  obtained  by  least  squares  fit  to  the  computed  values  is  about  —0.4.  Note  that 
the  effective  plastic  strain-rate  at  the  crack-tip  is  hnite,  and  the  plot  in  Figure  4.9  is  for 
r/a  >  10“^'^,  where  a  is  the  current  crack  length.  Results  plotted  in  Figure  4.9b  indicate 
that  without  cracking,  once  an  ASB  has  developed,  the  spatial  variation  of  the  effective 
plastic  strain-rate  ahead  of  the  ASB-tip  varies  very  slowly  and  its  value  at  the  ASB-tip 
equals  ~  1.31  x  10^/s.  It  does  not  exhibit  the  singularity  prevalent  in  the  effective  plastic 
strain-rate  ahead  of  the  crack-tip. 

We  have  plotted  in  Figure  4.10a-c  the  distribution  of  the  temperature  rise  in  W,  NiFe 
and  W2NiFe  plates  at  t  =  40,  170  and  70  /xs  respectively.  These  times  correspond  to  instants 
when  the  ASB/crack  has  propagated  ~  1  mm  from  the  point  of  initiation.  Temperatures 
at  the  ASB/crack-tips  at  these  times  are  listed  below.  Thus  the  temperature  rise  at  the 
crack-tip  is  nearly  14  K  higher  than  that  at  the  ASB-tip.  For  the  case  of  no  crack  opening. 
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Table  4.3:  Temperature  rise  at  ASB/crack-tip 


Material 

ASB-tip 

(K) 

Crack-tip 

(K) 

W 

512 

526 

NiFe 

547 

561 

W2NiFe 

546 

561 

temp:  150.0  230.0 


Figure  4.10:  Temperature  distribution  in  (a)  W,  (b)  NiFe  and  (c)  W2NiFe  plates  at  t  = 
40,  170  and  70  /rs,  respectively.  Plots  on  the  left  (right)  are  without  (with) 
crack  opening 


the  maximum  temperature  occurs  at  the  specimen  centroid  and  not  at  the  ASB-tip.  This  is 
because,  with  the  passage  of  time,  the  difference  between  the  effective  plastic  strain  at  the 
specimen  centroid  from  where  an  ASB  hrst  originated  and  that  at  the  ASB-tip  continues 
to  increase;  e.g.  see  Figure  4.9b.  However,  when  a  crack  is  formed,  then  the  maximum 
temperature  occurs  at  the  crack-tip.  Values  listed  in  Table  4.3  depend  upon  the  ASB/crack 
initiation  criterion. 

Relative  sliding  between  surfaces 

The  relative  sliding  between  the  upper  and  the  lower  parts  of  the  body  is  exhibited  in  Figure 
4.11  that  shows  the  deformed  shape  of  a  small  region.  It  is  clear  that  the  part  of  the  body 
occupying  the  region  X2  >  0  in  the  reference  conhguration  has  moved  to  the  right  relative 
to  that  in  the  domain  Xi  <  0.  Recalling  that  the  FE  mesh  in  the  reference  conhguration 
has  rectangular  elements  with  sides  parallel  to  the  horizontal  and  the  vertical  axes,  the 
two  elements  in  the  vertical  direction  next  to  the  sliding/crack  surface  have  been  sheared 
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considerably  more  than  those  further  away  from  the  crack  surface.  The  slight  asymmetry 
in  the  results  about  the  horizontal  centroidal  axis  is  attributed  to  numerical  oscillations  in 
the  computed  results.  The  smoothing  of  oscillations  located  the  crack  surface  at  X2  =  0.001 
rather  than  at  X2  =  0.0.  At  some  points,  the  crack  surface  X2  =  0+  has  moved  over  the 
crack  surface  X2  =  0“  by  as  much  as  four  element  lengths,  or  ~  18/im. 


0.6  0.62  0.64  0.66 

X 


Figure  4.11:  A  section  of  the  W  plate  deformed  in  plane  strain  shear  showing  relative 
sliding  between  the  cracked  surfaces 


Porosity  Evolution 

The  hydrostatic  pressure  in  the  material  within  the  ASB,  ahead  of  it  and  directly  in  front  of 
the  crack-tip  was  found  to  be  compressive.  Thus  no  new  voids  nucleated  in  the  material  even 
though  it  had  enormous  plastic  deformations.  The  total  porosity  evolved  in  the  intensely 
deformed  region  is  rather  minuscule. 

4.3.5  Dependence  of  Results  on  FE  Mesh 

In  Section  4.3.1  we  showed  that  an  element  size  of  100  /rm  x  100  /rm  for  a  FG  linear  elastic 
bar  gives  numerical  results  close  to  the  analytical  solution  of  the  wave  propagation  problem. 
As  shown  in  Table  3.1,  the  uniform  FE  meshes  of  element  size  125  /rm  x  125  fim  and  42  fim 
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X  42  /iin  gave  ASB  initiation  times  differing  by  2.1%.  Thus,  FE  meshes  used  herein  should 
give  reasonably  accurate  results.  The  CPU  time  needed  to  compute  speeds  of  propagation  of 
an  ASB  and  a  crack  due  to  ductile  failure  in  a  FG  body  is  close  to  200  hours.  Thus,  results 
could  not  be  computed  with  a  FE  mesh  much  hner  than  the  one  used  here. 

One  way  to  eliminate  the  effect  of  FE  mesh  on  computed  results  is  to  use  a  strain-rate 
gradient  dependent  viscoplasticity  theory  such  as  that  employed  in  [97-99].  It  introduces 
material  characteristic  lengths  that  can  not  be  easily  estimated. 

4.4  Conclusions 

We  have  analyzed  the  initiation  and  propagation  of  brittle  and  ductile  fractures  in  homo¬ 
geneous  and  functionally  graded  plates  deformed  either  in  plane  strain  tension  at  nominal 
strain-rates  of  200  and  2,000/s  or  in  plane  strain  shear  at  a  nominal  strain-rate  of  5,000/s. 
It  is  found  that  for  the  crack  propagating  due  to  brittle  failure  at  the  nominal  strain-rate  of 
200/s,  the  driving  force  is  a  non-decreasing  function  of  the  crack  extension.  At  the  higher 
strain-rate  of  2,000/s  the  axial  pulling  force  increases  as  the  crack  extends  till  the  crack 
length  equals  approximately  one-half  the  plate  width,  and  it  subsequently  decreases.  When 
the  crack  length  in  a  NiFe  plate  deformed  at  an  average  strain-rate  of  2,000/s  equals  about 
one-half  the  plate  width,  a  large  region  of  the  material  ahead  of  the  crack-tip  fails  instan¬ 
taneously  signifying  shattering  of  the  plate.  It  is  found  that  the  NiFe  plate  shatters  only  if 
the  nominal  axial  strain-rate  exceeds  1,130/s,  and  the  crack  elongation  prior  to  shattering 
varies  with  the  nominal  axial  strain-rate.  For  the  hnite  size  plate  studied  here,  the  maximum 
computed  crack  speed  in  a  FG  plate  pulled  at  a  nominal  axial  strain-rate  of  2,000/s  is  almost 
equal  to  the  Rayleigh  wave  speed.  For  an  FG  plate,  the  Rayleigh  wave  speed  varies  with  the 
position. 

We  have  also  analyzed  ductile  failure  in  a  plate  deformed  in  plane  strain  simple  shear. 
It  is  found  that  in  a  W  plate,  an  adiabatic  shear  band  (ASB)  and  a  crack  propagate  at 
virtually  identical  speeds.  However,  in  a  NiFe  plate  the  crack  speed  is  higher  than  the  ASB 
speed  in  the  beginning  but  the  two  are  nearly  equal  after  they  have  propagated  for  15  /rs.  In 
W  and  NiFe  plates,  the  crack  speed  increases  from  ~  0.05  km/s  at  the  instant  of  initiation 
to  ~  1.8  km/s  when  the  crack-tip  has  approached  the  edge  of  the  plate.  However,  in  a  FG 
plate  with  the  crack  propagating  from  W  rich  to  W  poor  regions,  the  crack  propagation 
speed  is  nearly  steady  for  a  certain  interval,  and  its  maximum  value  is  about  1/lOth  of  that 
in  a  homogeneous  W  or  NiFe  plate.  The  tangential  force  applied  at  the  top  and  the  bottom 
surfaces  decreases  as  the  ASB/crack  extends,  and  is  lower  when  a  crack  is  allowed  to  open 
than  that  without  the  opening  of  a  crack.  The  effective  plastic  strain-rate  is  hnite  at  the 
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crack-tip  but  is  singular  ahead  of  the  propagating  crack-tip;  the  order  of  singularity  is  ~  0.4. 
However,  the  effective  plastic  strain-rate  ahead  of  a  propagating  ASB  varies  rather  gradually 
and  equals  ~  1.3  x  10^/s  in  a  W  plate  deformed  at  a  nominal  strain-rate  of  5,000/s.  The 
temperature  at  a  crack-tip  is  ~  14  K  higher  than  that  at  an  ASB-tip. 


68 


Chapter  5 


Comparison  of  Results  from  Meso- 
and  Macro-level  Computations  of 
Adiabatic  Shear  Bands  and 
Multi-scale  Modeling 

5.1  Problem  Statement 

Recall  that  in  Chapter  2,  we  studied  plane  strain  deformations  of  a  RVE  and  found  val¬ 
ues  of  material  parameters  for  the  homogeneous  medium  so  that  the  responses  of  the  RVE 
and  the  homogenized  medium  in  thermal  and  mechanical  deformations  were  virtually  the 
same.  Naturally,  the  following  question  arises:  will  an  ASB  initiate  at  the  same  time  when 
deformations  of  the  individual  constituents  are  modeled  as  compared  to  the  case  when  de¬ 
formations  of  the  equivalent  homogenized  medium  are  scrutinized?  In  the  former  case,  the 
two  materials  are  presumed  to  be  perfectly  bonded  together  at  their  common  interfaces. 
Because  of  differences  in  their  acoustic  and  thermal  impedances,  mechanical  and  thermal 
waves  are  reflected  and  refracted  from  the  interfaces  between  the  two  materials.  Depending 
upon  the  number  of  particulates  in  the  composite,  these  waves  are  reflected  and  refracted 
numerous  times.  However,  in  the  homogenized  material  waves  are  only  reflected  from  the 
outer  bounding  surfaces.  The  interaction  among  these  waves  and  also  their  interaction  with 
the  incident  waves  result  in  an  inhomogeneous  deformation  of  the  body.  Below,  we  delineate 
the  influence  of  the  particulate  size  on  the  initiation  and  propagation  of  an  ASB,  and  then 
compare  the  formation  of  an  ASB  with  that  in  the  equivalent  homogenized  medium. 


69 


5.2  ASBs  in  a  Particulate  Composite 


As  in  Chapters  3  and  4,  plane  strain  coupled  thermomechanical  deformations  of  the  material 
in  the  hrst  quadrant  are  analyzed.  The  0.5  mm  x  0.5  mm  region  is  divided  into  100  x  100 
uniform  4-node  quadrilateral  elements.  Keeping  the  volume  fraction  of  W  particulates  hxed 
at  50%,  random  distributions  of  particulates  of  40,  50,  60  and  75  /xm  in  diameter,  as  shown 
in  Figure  5.1,  are  considered.  Boundary  conditions  arising  from  the  presumed  symmetry  of 
deformations  are  applied  on  the  bottom  and  the  left  edges,  and  the  other  two  edges  are  taken 
to  be  thermally  insulated.  The  right  edge  is  traction  free  and  the  prescribed  normal  velocity 
on  the  top  smooth  surface  increases  linearly  from  0  to  2.5  m/s  in  1  /xs  and  then  stays  steady. 
Thus  the  nominal  axial  strain-rate  induced  in  the  body  is  5000/s.  Figure  5.1  also  exhibits 
time-histories  of  the  axial  loads  for  the  four  cases.  As  for  a  homogeneous  body,  with  the 
continuous  application  of  the  axial  velocity,  the  axial  load  hrst  increases,  reaches  a  plateau 
and  then  drops  rapidly.  Whereas  the  time  history  of  the  axial  load  exhibits  oscillations  for  a 
homogeneous  body  (cf.  Figure  3.3),  it  varies  rather  smoothly  for  the  particulate  composite. 
It  suggests  that  the  overall  response  to  dynamic  loads  of  a  body  containing  a  reasonably 
large  number  of  particulates  is  smoother  than  that  of  a  homogeneous  body.  Even  though 
the  peak  axial  load  varies  slightly  with  the  particulate  size,  it  occurs  approximately  at  an 
axial  strain  of  0.08  for  the  four  particulate  sizes  considered. 

In  order  to  delineate  the  initiation  and  formation  of  an  ASB  in  a  particulate  composite, 
we  have  plotted  in  Figure  5.2  contours  of  the  effective  plastic  strain  at  hve  different  times, 
for  the  case  with  40  /rm  diameter  particulates.  As  the  top  surface  is  pulled  axially  with 
the  prescribed  velocity,  the  softer  NiFe  matrix  undergoes  more  severe  plastic  deformations 
as  compared  to  those  of  W  particulates.  In  Figure  5.2a,  there  are  several  narrow  regions 
of  intense  plastic  deformation.  During  deformations  of  the  body  from  t  =  30  /xs  to  t  =  40 
/is,  the  effective  plastic  strain  has  grown  in  many  of  these  regions.  However,  it  is  hard  to 
delineate  where  an  ASB  will  eventually  develop.  There  seems  to  be  at  least  two  regions  of 
large  plastic  deformation  in  Figures  5.2c  and  5. 2d.  Even  at  an  average  axial  strain  of  0.35 
there  are  many  regions  where  the  effective  plastic  strain  is  nearly  100%.  We  now  attempt  to 
delineate  which  ones,  if  any,  of  these  regions  correspond  to  ASBs. 

Figure  5.3  shows  the  axial  stress  ayy  and  the  effective  plastic  strain  as  a  function  of 
time  in  the  40  /rm  particulate  diameter  composite  at  two  points:  a  point  inside  W  initially 
located  at  {Xi,X2)  =  (0.06,0.4)  and  a  nearby  point  located  inside  the  NiFe  matrix  initially 
at  (0.085,0.4).  Note  that  the  effective  plastic  strain  is  much  higher  in  the  NiFe  matrix  than 
in  the  W  particulate  and  the  axial  stress  is  higher  in  the  W  than  in  the  NiFe.  Because  of 
the  assumption  of  perfect  bonding  between  the  two  materials,  displacements  and  surface 
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(e) 

Figure  5.1:  W  (particulate) /NiFe  (matrix)  RVEs  with  particulate  diameter  (a)  40  /xm, 
(b)  50  ^m,  (c)  60  /xm,  and  (d)  75  ^m;  (e)  axial  load  versus  time  for  RVEs 
(a)-(d) 


tractions  at  common  interfaces  between  them  are  continnons.  Even  thongh  the  qnasistatic 
yield  stress  of  W  is  nearly  three  times  that  of  NiFe,  for  t  >  10  /xs  the  effective  stress  in  W 
is  only  a  little  higher  than  that  in  NiFe  becanse  of  the  large  differences  in  their  strain-  and 
strain-rate  hardening  characteristics.  However,  results  plotted  in  Figure  5.3b  reveal  that  the 
effective  plastic  strain  in  W  is  less  than  one-fourth  that  in  NiFe.  Since  the  heat  capacity 
pc  of  W  is  approximately  24%  higher  than  that  of  NiFe,  the  temperature  rise  is  W  is  only 
slightly  less  than  that  of  NiFe;  this  is  depicted  in  Figure  5.3c.  Thus  the  thermal  energy 
density,  pc0,  in  W  is  nearly  20%  less  than  that  in  NiFe.  The  two  materials  are  softened 
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thermally  at  different  rates  due  to  differences  in  their  temperatures  and  thermal  softening 
parameters.  Thus  effective  plastic  strain  rates  are  different  in  the  two  materials.  Results 
plotted  in  Figure  5.3b  suggest  that  the  effective  plastic  strain  rate  at  the  point  in  the  NiFe 
matrix  is  considerably  more  than  that  in  the  W  particulate.  At  t  ~  70  /rs  or  the  average 
axial  strain  of  about  0.35,  the  effective  plastic  strain  rate  at  the  point  in  NiFe  decreases 
noticeably. 

For  t  =  80.1  /is  and  t  =  90.1  /xs,  Figure  5.4  depicts  points  where  the  ASB  initiation 
criterion  adopted  in  Section  3.2.2  has  been  satished  and  also  contours  of  the  effective  plastic 
strain.  Points  where  the  ASB  criterion  has  been  satished  do  not  lie  on  a  continuous  curve 
and  hence  do  not  form  a  discernible  ASB,  despite  the  high  levels  of  deformation  exhibited 
in  the  corresponding  contours  of  the  effective  plastic  strain.  It  is  thus  clear  that  the  ASB 
initiation  criterion  in  a  particulate  composite  is  different  from  that  in  a  body  with  either 
uniform  or  continuous  variation  of  material  properties. 

After  a  close  scrutiny  of  the  spatial  and  the  temporal  variations  of  different  measures 
of  deformation,  it  was  found  that  at  late  stages  of  deformation  the  spatial  variations  of 
velocity  and  the  rate  of  change  of  temperature  are  similar.  At  points  where  the  velocity 
varies  sharply,  the  rate  of  increase  of  temperature  is  also  very  high.  In  locally  adiabatic 
deformations,  the  rate  of  increase  of  temperature  equals  the  energy  dissipation  rate/mass 
divided  by  the  specihc  heat.  In  general,  the  heat  conduction,  the  energy  dissipation  rate  and 
the  heat  capacity  influence  the  rate  of  temperature  rise  at  a  point. 

For  a  particulate  composite  with  50%  volume  fraction  of  40  /rm  diameter  W  particulates 
in  NiFe  matrix.  Figure  5.5  evinces  contour  plots  of  the  axial  velocity  and  the  rate  of  increase 
of  temperature  at  an  axial  strain  of  35%.  It  is  clear  that  the  axial  velocity  varies  sharply 
across  regions  where  the  rate  of  increase  of  temperature  is  also  extremely  large.  These  results 
and  those  plotted  in  Figure  5.3b  suggest  that  an  ASB  initiates  at  a  point  in  the  NiFe  matrix 
at  an  average  axial  strain  of  ~0.35. 

One  would  like  to  replace  the  particulate  composite  with  a  homogeneous  material  with 
effective  properties  obtained  from  the  technique  described  in  Chapter  2.  However,  the  par¬ 
ticulate  composite  possesses  many  particulate/matrix  interfaces  that  act  as  local  defects, 
producing  ASBs  earlier  than  would  be  expected  in  a  homogeneous  material.  We  propose 
to  simulate  the  cumulative  effect  of  these  defects  in  the  homogeneous  body  with  an  initial 
porosity  distribution  given  by  equation  5.1. 


fo{X) 


/ce  (1  -  #)  ,  0  <  r  < 
0,  r  >  H 


(5,1) 


Here  r  =  VX„X„,  is  the  distance  from  the  centroid  in  the  reference  conhguration  and 
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Figure  5.2:  Contours  of  effective  plas¬ 
tic  strain  for  an  RVE  with  W  particu¬ 
lates  of  diameter  40  /xm  at  (a)  30  /xs,  (b) 
40  /xs,  (c)  50  /xs,  (d)  60  /xs,  and  (e)  70  /xs 
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Figure  5.3:  Time  histories  of  (a)  axial 
stress,  (b)  effective  plastic  strain,  and  (c) 
temperature  rise  at  two  nearby  points, 
one  inside  a  particulate  and  the  other 
in  the  matrix  for  the  40  /xm  particulate 
diameter  composite 


H  =  0.5  mm  is  the  height  of  the  specimen.  The  porosity  fee  at  the  centroid  is  varied  such 
that  the  rate  of  increase  of  temperature  inside  the  ASB  in  the  particulate  composite  and  the 
homogeneous  material  are  similar.  Figure  5.6  shows  the  rate  of  increase  in  temperature  as  a 
function  of  time  for  several  points  inside  the  ASB  in  the  particulate  composites  and  the  rate 
of  increase  in  temperature  at  the  centroid  of  an  equivalent  homogeneous  body  for  different 
values  of  fee-  The  points  in  the  particulate  composites  were  selected  in  local  regions  of  intense 
plastic  deformation  that  also  have  a  large  rate  of  change  of  temperature.  Note  that  the  rate 
of  increase  of  temperature  for  the  homogeneous  body  increases  monotonically  with  time; 
that  for  the  particulate  composites  show  some  fluctuation  at  different  points.  Furthermore, 
at  any  time  t  several  points  in  a  particulate  composite  exhibit  heating  at  similar  rates.  This 
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Figure  5.4:  Points  satisfying  the  ASB  initiation  criterion  of  Section  3.2.2  and  contours 
of  effective  plastic  strain  in  a  W/NiFe  particulate  composite  with  50 
diameter  tungsten  particulates  at  (a)  80.1  fis  and  (b)  90.1  /rs 


behavior  is  not  seen  in  the  homogeneons  body,  where  the  centroid  is  heated  more  rapidly  than 
other  points  in  the  body.  These  severely  heated  points  in  a  particnlate  composite  may  be 
disconnected.  In  a  homogeneons  body,  the  ASB  originates  at  the  centroid  and  propagates 
at  45°  in  the  present  conhgnration.  ASB  formation  and  propagation  in  the  particulate 
composite  is  more  complex-ASBs  may  originate  either  simultaneously  or  sequentially  at 
several  points  and  then  connect,  forming  the  noticeable  shear  band. 

Dai  et  ah  [13]  studied  ASB  initiation  in  a  composite  body  comprised  of  SiC  particles  im¬ 
mersed  in  aluminum.  The  specimens  were  twisted  in  a  split  Hopkinson  bar  with  presumably 
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Figure  5.5:  Contours  of  (a)  axial  velocity  and  (b)  rate  of  increase  of  temperature  in  a 
particulate  composite  with  40  iim  diameter  particulates  at  an  axial  strain  of 
35% 


the  same  torsional  impulse.  They  found  that  the  maximum  shear  strain  induced  in  particu¬ 
late  composites  increased  with  a  decrease  in  the  particulate  diameter.  They  explained  this 
by  conjecturing  that  strain  gradients  increase  with  a  decrease  in  the  particulate  diameter 
and  this  provides  a  driving  force  for  the  ASB  initiation.  Batra  and  Kim  [100]  numerically 
simulated  displacement  controlled  simple  shearing  deformations  of  simple  and  dipolar  ma¬ 
terials  and  showed  that  the  consideration  of  strain-rate  gradients  and  the  corresponding 
higher-order  stresses  delays  the  initiation  of  an  ASB.  In  order  to  see  how  strain  gradients 
depend  upon  the  particulate  size,  we  have  plotted  in  Figure  5.7  the  spatial  variation  on 
the  line  X2  =  0.4  mm  of  the  effective  plastic  strain  for  the  four  particulate  sizes  considered 
when  the  average  axial  strain  in  each  composite  equals  0.25.  These  plots  conhrm  Dai  et 
al.’s  conjecture.  However,  one  can  not  conclude  from  our  simulations  that  ASBs  initiate 
sooner  with  a  decrease  in  the  particulate  diameter.  We  note  that  the  composite  studied  by 
Dai  et  ah  is  comprised  of  SiC  particles  disbursed  in  2024  aluminum  and  most  likely  the 
SiC  particles  are  hardly  deformed,  though  this  is  not  clear  from  the  micrographs  included 
in  [13].  Furthermore,  the  SiC  particles  considered  had  diameters  of  3.5  pm,  10  pm,  and  20 
pm  and  were  present  in  very  small  volume  fractions.  Zhou  [101]  has  analyzed  ASBs  in  plane 
strain  simple  shear  deformations  of  composites  comprised  of  75.4%  W  particulates  in  a  NiFe 
matrix  and  considered  particulate  diameters  of  10  pm,  5  pm,  and  2.5  pm.  A  sharp  notch 
provided  a  nucleation  site  for  an  ASB.  It  was  found  that  the  ASB  initiation  time  decreased 
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Figure  5.6:  Rate  of  increase  in  temperature  at  a  point  inside  the  ASB  in  50/50  W/NiFe 
particulate  composites  with  particulate  diameters  of  (a)  40  /um,  (b)  50  /xm, 
(c)  60  /xm,  and  (d)  75  /xm,  as  well  as  the  rate  of  increase  of  temperature  at  the 
centroid  of  equivalent  homogeneous  bodies  with  nonuniform  initial  porosity; 
coordinates  denote  the  position  of  the  point  in  the  particulate  composite  in 
the  reference  configuration. 


with  a  decrease  in  the  particulate  diameter  because  smaller  particles  facilitated  the  turning 
of  an  ASB  into  the  relatively  weak  matrix. 

In  order  to  assess  the  effect  of  the  particulate  arrangement,  we  have  plotted  in  Figure  5.8 
contour  plots  of  the  axial  velocity  for  three  different  random  distributions  of  50  /xm  diameter 
particulates.  It  is  clear  that  the  particulate  placements  strongly  influence  when  and  where 
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an  ASB  initiates.  For  example  in  Figure  5.8a,  an  ASB  seems  to  have  originated  from  a  point 
on  the  vertical  centroidal  axis  and  propagated  at  ±45°,  but  in  Figure  5.8c,  the  origin  of  the 
ASB  is  at  a  point  on  the  free  surface  and  it  propagated  only  towards  the  centroid  of  the 
cross-section.  Furthermore,  the  times  of  initiation  equal  ~90  /xs  and  60  /xs  for  particulate 
distributions  of  Figures  5.8a  and  5.8c,  respectively.  For  the  particulate  distribution  of  Figure 
5.8b,  there  is  no  narrow  region  across  which  the  axial  velocity  increases  sharply.  The  plots  of 
the  rate  of  temperature  increase  depicted  in  Figure  5.9  support  the  aforementioned  remarks. 
The  variation  of  the  effective  plastic  strain  on  the  line  X2  =  0.4  mm  is  exhibited  in  Figure 
5.10.  It  is  apparent  that  strain  gradients  are  higher  for  the  particulate  arrangement  of  Figure 
5.10b  than  those  for  the  other  two  arrangements  analyzed.  The  time  histories  of  the  rate  of 
temperature  rise  evinced  in  Figure  5.11  suggest  that  an  ASB  initiates  at  t  ~  60  /xs.  There 
is  less  variation  in  the  ASB  initiation  time  than  in  the  contour  plots  of  the  axial  velocity 
and  the  rate  of  temperature  increase.  An  initial  affine  porosity  variation  with  approximately 
1.5%  porosity  at  the  specimen  centroid  in  the  equivalent  homogenized  body  will  result  in 
the  ASB  initiating  at  t  ~  60  /is. 

We  note  that  the  time  histories  of  the  axial  load  plotted  in  Figure  5.1  indicate  that  the 
drop  in  the  axial  load  to  80%  of  its  peak  value  occurs  at  64  /xs  <  t  <  68  /rs  for  the  four 
particulate  diameters.  Whether  or  not  it  can  be  used  as  a  criterion  for  the  ASB  initiation 
remains  to  be  investigated. 

The  particulate  diameters  considered  here  are  typical  of  those  employed  by  Dick  et  ah 
[102]  in  their  experiments  on  tungsten  heavy  alloys. 


78 


effective  plastic  strain  effective  plastic  strain 


Figure  5.7:  At  an  average  axial  strain  of  0.25,  variation  of  the  effective  plastic  strain 
along  X2  =  0.4  mm  for  particulate  composites  with  particulate  diameters  of 
(a)  40  ^m,  (b)  50  fiui,  (c)  60  ^um,  and  (d)  75  fiui 
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Figure  5.8:  Contours  of  the  axial  veloc¬ 
ity  for  three  different  random  arrange¬ 
ments  of  50  ^m  diameter  particulates: 
(a)  arrangement  1  at  90  ^s,  (b)  arrange¬ 
ment  2  at  70  fj,s,  and  (c)  arrangement  3 
at  60  jjs. 
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Figure  5.9:  Contours  of  the  rate  of 
change  of  temperature  for  three  different 
random  arrangements  of  50  /xm  diame¬ 
ter  particulates:  (a)  arrangement  1  at  90 
jjs,  (b)  arrangement  2  at  70  /xs,  and  (c) 
arrangement  3  at  60  /xs. 
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Figure  5.11:  Rate  of  increase  in  tem¬ 
perature  at  a  point  inside  the  ASB  in 
50/50  W/NiFe  particulate  composites 
with  particulate  diameters  of  50  ;um 
for  (a)  arrangement  1,  (b)  arrangement 
2,  and  (c)  arrangement  3,  as  well  as 
the  rate  of  increase  of  temperature  at 
the  centroid  of  equivalent  homogeneous 
bodies  with  nonuniform  initial  porosity; 
coordinates  denote  the  position  of  the 
point  in  the  particulate  composite  in  the 
reference  configuration. 
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5.3  Multiscale  Modeling 

We  first  analyze  the  initiation  and  development  of  an  ASB  in  a  prismatic  body  of  10  x  10 
mm  square  cross-section  deformed  in  plane  strain  tension.  The  body  is  made  of  a  homoge¬ 
nized  material  whose  response  is  equivalent  to  that  of  a  50%  volume  fraction  W  particulate 
composite  with  50  /rm  diameter  W  particulates  in  the  NiFe  matrix;  material  parameters  of 
the  homogenized  body  are  ascertained  by  using  the  procedure  outlined  in  Chapter  2.  As  in 
the  previous  studies  discussed  in  Chapter  3,  symmetries  about  the  two  centroidal  axes  are 
exploited  and  coupled  thermomechanical  deformations  of  the  material  in  the  hrst  quadrant 
are  analyzed.  The  5  mm  x  5  mm  region  is  divided  into  a  uniform  100  x  100  FE  mesh  of 
4-node  quadrilateral  elements.  The  initial  porosity  distribution  is  given  by  equation  (5.1) 
with  fee  =  0.025,  and  it  is  deformed  at  a  nominal  strain-rate  of  5,000/s  till  plastic  deforma¬ 
tions  have  localized  into  a  narrow  region.  The  time  histories  of  13  unknowns  at  nodes  on  the 
surfaces  Xi  =  1  mm  and  X2  =  1  mm  in  the  undeformed  conhguration  are  recorded  for  later 
use.  Subsequently,  detailed  deformations  of  the  material  in  the  1  mm  x  1  mm  region  with 
two  edges  along  the  centroidal  axes  of  the  original  cross-section  are  analyzed.  As  shown  in 
Figure  5.12,  the  0.5  mm  x  0.5  mm  region  has  50  /rm  diameter  W  particulates  randomly  dis¬ 
tributed  in  the  NiFe  matrix  and  the  remaining  region  is  the  homogenized  material  equivalent 
to  the  particulate  composite.  On  the  surfaces  Xi  =  1  mm  and  X2  =  1  mm  the  time  history 
of  surface  tractions  (normal  and  tangential  stresses)  and  of  temperature  as  computed  from 
the  previous  analysis  is  applied  with  the  remaining  two  surfaces  having  the  same  boundary 
conditions  as  in  the  original  analysis  of  the  5  mm  x  5  mm  region.  Numerical  experiments 
with  the  velocity  and  the  temperature  prescribed  on  the  surfaces  Xi  =  1  mm  and  X2  =  1  mm 
failed  because  of  severe  differences  in  deformations  of  particles  on  the  common  interface  be¬ 
tween  the  homogenized  medium  and  the  particulate  composite.  However,  these  experiments 
were  successful  for  a  homogenized  body  and  the  two  approaches  yielded  identical  results. 
With  tractions  assigned  on  the  surface  Xi  =  1  mm  and  X2  =  1  mm,  deformations  of  these 
and  other  particles  adopt  readily  to  the  constraints  imposed  by  the  stiffer  W  particulates.  It 
should  be  added  that,  in  the  reference  configuration,  nodes  on  the  top  and  the  right  surfaces 
of  the  1  mm  x  1  mm  mesh  coincide  with  those  on  the  surfaces  Xi  =  1  mm  and  X2  =  1  mm 
in  the  5  mm  x  5  mm  mesh,  respectively. 

5.3.1  Verification  of  the  Methodology 

The  aforementioned  approach  was  used  to  analyze  deformations  of  a  homogeneous  body. 
Figure  5.13  compares  the  time  histories  of  the  temperature  and  rate  of  increase  of  tempera¬ 
ture  of  a  material  particle  located  at  the  specimen  centroid  obtained  by  the  two  procedures. 
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Figure  5.12:  FE  mesh  for  the  sub-section  in  the  multiscale  analysis 


and  Figure  5.14  depicts  contours  of  the  effective  plastic  strain  at  an  average  axial  strain  of 
0.29  in  the  macro-scale  analysis.  It  is  clear  that  the  average  axial  strain  in  the  1  mm  x  1  mm 
subsection  is  more  than  0.29.  The  rate  of  increase  of  the  temperature  is  much  higher  in  the 
micro-macro  analysis  than  in  the  pure  macro  analysis.  The  combined  analysis  can  be  viewed 
as  rehning  the  FE  mesh  and  as  shown  by  Batra  and  Ko  [21]  the  effective  plastic  strain  rate 
at  the  ASB  center  increases  when  FE  mesh  is  rehned.  The  ASB  initiates,  as  signihed  by 
the  jump  in  the  rate  of  temperature  increase,  slightly  earlier  for  the  macro-micro  analysis  as 
compared  to  that  for  the  macro  analysis;  the  difference  between  the  two  times  is  ~7%.  The 
intensely  deformed  region  is  signihcantly  narrower  for  the  macro-micro  analysis  than  that 
for  the  macro  analysis.  In  Eigure  5.14a,b  the  right  edge  is  not  straight  because  it  is  an  inner 
surface  in  the  5  mm  x  5  mm  cross-section.  For  the  same  reason,  the  average  axial  strain  in 
this  1  mm  x  1  mm  region  need  not  equal  the  average  axial  strain  of  0.29  in  the  5  mm  x  5 
mm  region.  This  micro-macro  analysis  could  be  repeated  several  times  to  fully  resolve  an 
ASB. 
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(a)  (b) 

Figure  5.13:  (a)  Temperature  and  (b)  rate  of  increase  of  temperature  versus  time  at  the 
centroid  of  a  homogenized  (50%  W,  50%  NiFe)  body  for  the  macro  and  the 
macro-micro  analyses 
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Figure  5.14:  Contours  of  effective  plastic  strain  for  a  homogenized  50/50  W/NiFe  sample 
at  an  average  axial  strain  of  0.29  for  Xi  <  1  mm  and  X2  <  1  mm  for  the 
(a)  macro-scale  analysis  and  (b)  macro-micro  analysis 


5.3.2  Results  for  the  Particulate  Composite 


The  above  described  procedure  was  applied  to  four  different  distributions  of  50  /im  diameter 
W  particulates  disbursed  in  the  0.5  mm  x  0.5  mm  region.  Figure  5.15  depicts  contours  of 
the  effective  plastic  strain  and  Figure  5.16  contours  of  the  rate  of  increase  of  temperature 
for  the  four  distributions.  It  is  evident  that  an  ASB  initiates  from  a  different  point  in  each 
case.  Furthermore,  intense  plastic  deformations  occur  in  different  regions  in  the  four  cases, 
and  only  one  of  these  agrees  with  the  general  shape  of  the  intensely  deformed  region  in  the 
equivalent  homogeneous  body;  c.f.  Figure  5.9.  The  ASB  center  can  be  located  from  a  visual 
inspection  of  the  contour  plots  of  Figures  5.15  and  5.16;  its  location  does  not  coincide  with 
the  origin  or  the  centroid  of  the  cross-section,  where  the  initial  porosity  is  maximum.  The 
inhomogeneities  induced  by  the  presence  of  W  particulates  strongly  affect  where  and  when 
an  ASB  initiates.  For  the  four  distributions  of  W  particulates,  we  have  plotted  in  Figure 
5.17  time  histories  of  the  rate  of  temperature  rise  and  the  effective  plastic  strain  at  the  ASB 
center.  The  ordered  and  the  three  random  distributions  of  particulates  give  ASB  initiation 
times  between  22  and  24  /rs  which  is  considerably  sooner  than  the  nearly  60  /xs  deduced  from 
the  meso-scale  results  plotted  in  Figure  5.6b.  The  large  difference  between  the  two  values 
can  not  be  attributed  to  the  considerably  hner  mesh  employed  in  the  micro  analysis;  the 
reasons  are  a  subject  for  future  investigations.  The  deformed  conhguration  plotted  in  Figure 
5.15  suggest  that  the  nominal  axial  strain  for  the  three  random  distributions  is  essentially 
the  same  but  is  slightly  less  than  that  for  the  ordered  distribution  of  particulates. 
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Figure  5.15:  Contours  of  the  effective  plastic  strain  at  t  =  25  /rs  for  (a)  an  ordered 
particulate  arrangement  and  (b)-(d)  three  random  particulate  arrangements 
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Figure  5.16:  Contours  of  the  rate  of  increase  in  temperature  at  t  =  25  /rs  for  (a)  an 
ordered  particulate  arrangement  and  (b)-(d)  three  random  particulate  ar¬ 
rangements 
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(a) 


(b) 


Figure  5.17:  (a)  Effective  plastic  strain  and  (b)  rate  of  temperature  increase  at  a  point 
inside  the  ASB  as  functions  of  time  for  four  particulate  arrangements 
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5.4  Conclusions 


We  have  studied  the  initiation  and  propagation  of  ASBs  in  particulate  composites,  using 
both  plane  strain  tension  tests  at  the  meso-scale  and  a  multi-scale  analysis  that  bridges 
the  length  scales  between  meso-  and  macro-level  analyses  of  the  heterogeneous  composite.  It 
was  found  that  the  initiation  and  propagation  of  ASBs  in  particulate  composites  is  decidedly 
different  than  that  in  the  equivalent  homogenized  material  (with  properties  derived  from  the 
technique  presented  in  Chapter  2);  particulate/matrix  interfaces  act  as  defects,  promoting 
the  formation  of  ASBs.  We  have  modeled  these  defects  in  the  homogenized  material  by 
assuming  that  the  initial  porosity  varies  affinely  with  the  maximum  value,  fee,  at  the  centroid 
of  the  cross-section.  Values  of  fee  between  0.0075  and  0.015  give  ASB  initiation  times  for 
the  homogenized  body  that  are  close  to  those  for  the  particulate  composite. 

The  multiscale  analyses  of  the  particulate  composite  revealed  that  ASB  initiation  time 
is  nearly  one-half  of  that  for  the  meso-scale  analysis  of  a  composite  with  the  same  volume 
fraction  and  distribution  of  particulates.  Reasons  for  this  difference  are  unclear  and  need  to 
be  investigated. 

The  ASB  initiation  criterion  for  heterogeneous  solids  is  totally  different  from  that  for 
homogeneous  bodies. 
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Chapter  6 
Contributions 


A  significant  contribution  of  the  present  work  is  the  examination  of  a  particulate  composite 
under  a  wide  variety  of  loading  conditions  at  multiple  length  scales.  A  numerical  homog¬ 
enization  technique  for  thermoelastoviscoplastic  particulate  composites  was  developed,  en¬ 
abling  one  to  calculate  effective  viscoplastic  parameters  that  are  not  easily  obtained  through 
traditional  analytical  methods.  Adiabatic  shear  bands  were  studied  for  the  hrst  time  in  func¬ 
tionally  graded  materials,  demonstrating  that  FGMs  can  produce  quantitatively  and  quali¬ 
tatively  different  results  from  homogeneous  materials.  Furthermore,  the  brittle  and  ductile 
failure  of  FGMs  was  studied;  for  the  hrst  time,  a  sustained  ductile  crack  was  simulated 
and  analyzed.  The  initiation  and  propagation  of  adiabatic  shear  bands  at  the  microstruc- 
tural  level  was  also  studied-it  was  shown  that  prior  ASB  criteria  were  not  applicable  to  a 
multi-constituent  particulate  composite.  These  contributions  are  detailed  below. 

We  have  used  the  hnite  element  technique  to  study  plane  strain  deformations  of  repre¬ 
sentative  volume  elements  of  tungsten/nickel-iron  particulate  composites.  Through  a  con¬ 
structed  set  of  tests,  effective  material  parameters  such  as  Young’s  modulus,  Poisson’s  ratio, 
strain-  and  strain-rate  hardening  parameters,  thermal  softening  parameters,  thermal  conduc¬ 
tivity  and  the  coefficient  of  thermal  expansion  are  determined  for  each  given  microstructure. 
It  was  found  that  particulate  size  and  arrangement  played  little  role  in  effective  material 
properties  as  long  as  the  arrangement  is  statistically  isotropic.  The  effect  of  the  volume 
fraction  of  constituents  was  also  examined-effective  material  parameters  for  volume  frac¬ 
tions  between  50  and  80%  of  tungsten  particulate  were  obtained  for  tungsten/nickel-iron 
(frequently  referred  to  as  tungsten-heavy  alloy).  The  signihcant  contribution  is  the  homog¬ 
enization  technique,  however-it  can  easily  be  adapted  to  any  isotropic  composite,  giving 
effective  material  parameters  that  have  been  difficult  to  determine  previously. 

The  plane-strain  hnite  element  code  is  then  used  to  study  the  formation  of  adiabatic 
shear  bands  in  functionally  graded  materials.  FGMs  were  comprised  of  tungsten  and  nickel- 
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iron  (as  well  as  some  tungsten  and  copper /iron/steel  systems).  It  was  shown  that  ASBs 
tend  to  form  in  the  highly-ASB-susceptible  region  (in  this  case,  tungsten).  Furthermore,  it 
was  shown  that  ASB  initiation  occurred  sooner  in  systems  with  more  spatial  variation  of 
material  properties,  possibly  due  to  the  greater  degree  of  inhomogeneity. 

Brittle  and  ductile  failure  were  studied  in  tungsten  and  nickel-iron  as  well  as  in  W/NiFe 
FGMs.  Brittle  crack  speed  was  found  to  follow  the  Rayleigh  wave  speed  in  the  local  material 
given  a  sufficient  applied  nominal  strain-rate.  Ductile  failure  was  found  to  closely  match 
the  propagation  speed  of  an  ASB.  It  was  also  shown  that  ductile  failure  can  be  signihcantly 
delayed  and  the  ductile  crack  speed  altered  by  spatial  variation  in  material  properties.  To  the 
best  of  the  author’s  knowledge,  this  work  represents  the  hrst  time  the  sustained  propagation 
of  a  ductile  crack  has  been  simulated  numerically. 

A  multiscale  analysis  technique  was  used  to  study  the  initiation  and  propagation  of 
ASBs  in  particulate  composites.  Macro-scale  analyses  of  an  effective  homogeneous  material 
were  used  to  derive  boundary  conditions  for  a  micro-scale  analysis.  It  was  shown  that  ASB 
initiation  location  in  a  particulate  composite  is  dependent  on  the  microstructure  (such  as 
the  arrangement  of  particulates),  however,  the  initiation  time  is  not  signihcantly  affected. 
Furthermore,  it  was  shown  that  prior  ASB  initiation  criteria  are  not  applicable  to  particulate 
composites-load  transfer  between  the  constituents  makes  stress-based  criteria  unusable. 
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